Bug fix (wrong condition)

We should test the numbe of output arguments, not the number of input
arguments. This was bug was probably not affecting the outcomes since
the number of input arguments is always greater than 1.

doc/manual/source/conf.py

@@ -1,6 +1,6 @@
 # -*- coding: utf-8 -*-
 
-# Copyright © 2018-2023 Dynare Team
+# Copyright © 2018-2024 Dynare Team
 #
 # This file is part of Dynare.
 #
@@ -34,7 +34,7 @@ html_static_path = ['_static']
 master_doc = 'index'
 
 project = u'Dynare'
-copyright = u'1996–2023 Dynare Team'
+copyright = u'1996–2024 Dynare Team'
 author = u'Dynare Team'
 
 add_function_parentheses = False
@@ -71,7 +71,6 @@ latex_elements = {
                     warningBorderColor={RGB}{255,50,50},OuterLinkColor={RGB}{34,139,34}, \
                     InnerLinkColor={RGB}{51,51,255},TitleColor={RGB}{51,51,255}',
     'papersize': 'a4paper',
-    'preamble': r'\DeclareUnicodeCharacter{200B}{}', # Part of the workaround for #1707
 }
 
 latex_documents = [

doc/manual/source/dynare-misc-commands.rst

@@ -8,17 +8,17 @@
 Dynare misc commands
 ####################
 
-.. matcomm:: send_endogenous_variables_to_workspace
+.. matcomm:: send_endogenous_variables_to_workspace ;
 
     Puts the simulation results for the endogenous variables stored in ``oo_.endo_simul`` 
     into vectors with the same name as the respective variables into the base workspace.
 
-.. matcomm:: send_exogenous_variables_to_workspace
+.. matcomm:: send_exogenous_variables_to_workspace ;
 
     Puts the simulation results for the exogenous variables stored in ``oo_.exo_simul`` 
     into vectors with the same name as the respective variables into the base workspace.
 
-.. matcomm:: send_irfs_to_workspace
+.. matcomm:: send_irfs_to_workspace ;
 
     Puts the IRFs stored in ``oo_.irfs`` into vectors with the same name into the base workspace.
 

doc/manual/source/index.rst

@@ -11,7 +11,7 @@ Currently the development team of Dynare is composed of:
 * Willi Mutschler (University of Tübingen)
 * Johannes Pfeifer (University of the Bundeswehr Munich)
 * Marco Ratto (European Commission, Joint Research Centre - JRC)
-* Normann Rion (CY Cergy Paris Université and CEPREMAP)
+* Normann Rion (CEPREMAP)
 * Sébastien Villemot (CEPREMAP)
 
 The following people used to be members of the team:
@@ -26,7 +26,7 @@ The following people used to be members of the team:
 * Ferhat Mihoubi
 * George Perendia
 
-Copyright © 1996-2023, Dynare Team.
+Copyright © 1996-2024, Dynare Team.
 
 Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.3 or any later version published by the Free Software Foundation; with no Invariant Sections, no Front-Cover Texts, and no Back-Cover Texts.
 

doc/manual/source/the-model-file.rst

@@ -1259,7 +1259,8 @@ command, the list of transformed model equations using the
 ``write_latex_dynamic_model command``, and the list of static model
 equations using the ``write_latex_static_model`` command.
 
-.. command:: write_latex_original_model (OPTIONS);
+.. command:: write_latex_original_model ;
+             write_latex_original_model (OPTIONS);
 
     |br| This command creates two LaTeX files: one
     containing the model as defined in the model block and one
@@ -1334,7 +1335,8 @@ equations using the ``write_latex_static_model`` command.
         See :opt:`write_equation_tags`
 
 
-.. command:: write_latex_static_model (OPTIONS);
+.. command:: write_latex_static_model ;
+             write_latex_static_model (OPTIONS);
 
     |br| This command creates two LaTeX files: one
     containing the static model and one containing the LaTeX
@@ -1369,7 +1371,7 @@ equations using the ``write_latex_static_model`` command.
 
         See :opt:`write_equation_tags`.
 
-.. command:: write_latex_steady_state_model
+.. command:: write_latex_steady_state_model ;
 
     |br| This command creates two LaTeX files: one containing the steady
     state model and one containing the LaTeX document header
@@ -2945,8 +2947,8 @@ Finding the steady state with Dynare nonlinear solver
 
            ``5``
 
-                Newton algorithm with a sparse Gaussian elimination
+                Newton algorithm with a sparse Gaussian elimination (SPE)
-                (SPE) (requires ``bytecode`` option, see
+                solver at each iteration (requires ``bytecode`` option, see
                 :ref:`model-decl`).
 
            ``6``
@@ -2964,7 +2966,7 @@ Finding the steady state with Dynare nonlinear solver
            ``8``
 
                 Newton algorithm with a Stabilized Bi-Conjugate
-                Gradient (BICGSTAB) solver at each iteration (requires
+                Gradient (BiCGStab) solver at each iteration (requires
                 bytecode and/or block option, see :ref:`model-decl`).
 
            ``9``
@@ -3680,60 +3682,64 @@ speed-up on large models.
            ``0``
 
                Use a Newton algorithm with a direct sparse LU solver at each
-               iteration, applied on the stacked system of all the equations at
+               iteration, applied to the stacked system of all equations in all
-               every period (Default).
+               periods (Default).
 
            ``1``
 
                Use the Laffargue-Boucekkine-Juillard (LBJ) algorithm proposed
-               in *Juillard (1996)*. It is slower than ``stack_solve_algo=0``,
+               in *Juillard (1996)* on top of a LU solver. It is slower
-               but may be less memory consuming on big models. Note that if the
+               than ``stack_solve_algo=0``, but may be less memory consuming on
-               ``block`` option is used (see :ref:`model-decl`), a simple
+               big models. Note that if the ``block`` option is used (see
-               Newton algorithm with sparse matrices is used for blocks which
+               :ref:`model-decl`), a simple Newton algorithm with sparse
-               are purely backward or forward (of type ``SOLVE BACKWARD`` or
+               matrices, applied to the stacked system of all block equations
-               ``SOLVE FORWARD``, see :comm:`model_info`), since LBJ only makes
+               in all periods, is used for blocks which are purely backward or
-               sense on blocks with both leads and lags (of type ``SOLVE TWO
+               forward (of type ``SOLVE BACKWARD`` or ``SOLVE FORWARD``, see
-               BOUNDARIES``).
+               :comm:`model_info`), since LBJ only makes sense on blocks with
+               both leads and lags (of type ``SOLVE TWO BOUNDARIES``).
 
            ``2``
 
-               Use a Newton algorithm with a Generalized Minimal
+               Use a Newton algorithm with a Generalized Minimal Residual
-               Residual (GMRES) solver at each iteration (requires
+               (GMRES) solver at each iteration, applied on the stacked system
-               ``bytecode`` and/or ``block`` option, see
+               of all equations in all periods (requires ``bytecode`` and/or
-               :ref:`model-decl`)
+               ``block`` option, see :ref:`model-decl`)
 
            ``3``
 
-               Use a Newton algorithm with a Stabilized Bi-Conjugate
+               Use a Newton algorithm with a Stabilized Bi-Conjugate Gradient
-               Gradient (BICGSTAB) solver at each iteration (requires
+               (BiCGStab) solver at each iteration, applied on the stacked
-               ``bytecode`` and/or ``block`` option, see
+               system of all equations in all periods (requires ``bytecode``
-               :ref:`model-decl`).
+               and/or ``block`` option, see :ref:`model-decl`).
 
            ``4``
 
-               Use a Newton algorithm with an optimal path length at
+               Use a Newton algorithm with a direct sparse LU solver and an
-               each iteration (requires ``bytecode`` and/or ``block``
+               optimal path length at each iteration, applied on the stacked
-               option, see :ref:`model-decl`).
+               system of all equations in all periods (requires ``bytecode``
+               and/or ``block`` option, see :ref:`model-decl`).
 
            ``5``
 
-               Use a Newton algorithm with a sparse Gaussian
+               Use the Laffargue-Boucekkine-Juillard (LBJ) algorithm proposed
-               elimination (SPE) solver at each iteration (requires
+               in *Juillard (1996)* on top of a sparse Gaussian elimination
-               ``bytecode`` option, see :ref:`model-decl`).
+               (SPE) solver. The latter takes advantage of the similarity of
+               the Jacobian across periods when searching for the pivots
+               (requires ``bytecode`` option, see :ref:`model-decl`).
 
            ``6``
 
-               Synonymous for ``stack_solve_algo=1``. Kept for historical
+               Synonymous for ``stack_solve_algo=1``. Kept for backward
-               reasons.
+               compatibility.
 
            ``7``
 
-               Allows the user to solve the perfect foresight model
+               Allows the user to solve the perfect foresight model with the
-               with the solvers available through option
+               solvers available through option ``solve_algo``, applied on the
-               ``solve_algo`` (See :ref:`solve_algo <solvalg>` for a
+               stacked system of all equations in all periods (See
-               list of possible values, note that values 5, 6, 7 and
+               :ref:`solve_algo <solvalg>` for a list of possible values, note
-               8, which require ``bytecode`` and/or ``block`` options,
+               that values 5, 6, 7 and 8, which require ``bytecode`` and/or
-               are not allowed). For instance, the following
+               ``block`` options, are not allowed). For instance, the following
                commands::
 
                     perfect_foresight_setup(periods=400);
@@ -11739,7 +11745,7 @@ with ``discretionary_policy`` or for optimal simple rules with ``osr``
 
     With ``discretionary_policy``, the objective function must be quadratic.
 
-.. command:: evaluate_planner_objective;
+.. command:: evaluate_planner_objective ;
              evaluate_planner_objective (OPTIONS...);
 
     This command computes, displays, and stores the value of the
@@ -14293,7 +14299,7 @@ a trend target to which the endogenous variables may be attracted in the long-ru
    :math:`n\times 1` vector of parameters, :math:`A_i` (:math:`i=0,\ldots,p`)
    are :math:`n\times n` matrices of parameters, and :math:`A_0` is non
    singular square matrix. Vector :math:`\mathbf{c}` and matrices :math:`A_i`
-   (:math:`i=0,\ldots,p`) are set by Dynare by parsing the equations in the
+   (:math:`i=0,\ldots,p`) are set by parsing the equations in the
    ``model`` block. Then, Dynare builds a VAR(1)-companion form model for
    :math:`\mathcal{Y}_t = (1, Y_t, \ldots, Y_{t-p+1})'` as:
 
@@ -14504,7 +14510,7 @@ up to time :math:`t-\tau`, :math:`\mathcal{Y}_{\underline{t-\tau}}`) is:
 
 In a semi-structural model, variables appearing in :math:`t+h` (*e.g.*
 the expected output gap in a dynamic IS curve or expected inflation in a
-(New Keynesian) Phillips curve) will be replaced by the expectation implied by an auxiliary VAR
+New Keynesian Phillips curve) will be replaced by the expectation implied by an auxiliary VAR
 model. Another use case is for the computation of permanent
 incomes. Typically, consumption will depend on something like:
 
@@ -14512,13 +14518,13 @@ incomes. Typically, consumption will depend on something like:
 
       \sum_{h=0}^{\infty} \beta^h y_{t+h|t-\tau}
 
-Assuming that $0<\beta<1$ and knowing the limit of geometric series, the conditional expectation of this variable can be evaluated based on the same auxiliary model:
+Assuming that :math:`0<\beta<1` and knowing the limit of geometric series, the conditional expectation of this variable can be evaluated based on the same auxiliary model:
 
    .. math ::
 
       \mathbb E \left[\sum_{h=0}^{\infty} \beta^h y_{t+h}\Biggl| \mathcal{Y}_{\underline{t-\tau}}\right] = \alpha \mathcal{C}^\tau(I-\beta\mathcal{C})^{-1}\mathcal{Y}_{t-\tau}
 
-More generally, it is possible to consider finite discounted sums.
+Finite discounted sums can also be considered.
 
 .. command:: var_expectation_model (OPTIONS...);
 
@@ -15861,7 +15867,7 @@ Misc commands
     ``pdflatex`` and automatically tries to load all required
     packages. Requires the following LaTeX packages:
     ``breqn``, ``psfrag``, ``graphicx``, ``epstopdf``, ``longtable``,
-    ``booktabs``, ``caption``, ``float,`` ``amsmath``, ``amsfonts``,
+    ``booktabs``, ``caption``, ``float,`` ``amsmath``, ``amsfonts``, ``amssymb``, 
     and ``morefloats``.
 
 

doc/manual/utils/dynare_dom.py

@@ -1,6 +1,6 @@
 # -*- coding: utf-8 -*-
 
-# Copyright © 2018-2019 Dynare Team
+# Copyright © 2018-2024 Dynare Team
 #
 # This file is part of Dynare.
 #
@@ -80,9 +80,7 @@ class DynObject(ObjectDescription):
         signode += addnodes.desc_name(name, name)
 
         if self.has_arguments:
-            if not arglist:
+            if arglist:
-                signode += addnodes.desc_parameterlist()
-            else:
                 signode += addnodes.desc_addname(arglist,arglist)
         return fullname, prefix
 

license.txt

@@ -1,6 +1,6 @@
 Format: https://www.debian.org/doc/packaging-manuals/copyright-format/1.0/
 Upstream-Name: Dynare
-Upstream-Contact: Dynare Team, whose members in 2023 are:
+Upstream-Contact: Dynare Team, whose members in 2024 are:
                   - Stéphane Adjemian <stephane.adjemian@univ-lemans.fr>
                   - Michel Juillard <michel.juillard@mjui.fr>
                   - Frédéric Karamé <frederic.karame@univ-lemans.fr>
@@ -23,7 +23,7 @@ Upstream-Contact: Dynare Team, whose members in 2023 are:
 Source: https://www.dynare.org
 
 Files: *
-Copyright: 1996-2023 Dynare Team
+Copyright: 1996-2024 Dynare Team
 License: GPL-3+
 
 Files: matlab/+occbin/IVF_core.m
@@ -86,7 +86,7 @@ License: public-domain-aim
  Journal of Economic Dynamics and Control, 2010, vol. 34, issue 3,
  pages 472-489
 
-Files: matlab/optimization/bfgsi1.m matlab/csolve.m matlab/optimization/csminit1.m matlab/optimization/numgrad2.m
+Files: matlab/optimization/bfgsi1.m matlab/optimization/csolve.m matlab/optimization/csminit1.m matlab/optimization/numgrad2.m
  matlab/optimization/numgrad3.m matlab/optimization/numgrad3_.m matlab/optimization/numgrad5.m
  matlab/optimization/numgrad5_.m matlab/optimization/csminwel1.m matlab/+bvar/density.m
  matlab/+bvar/toolbox.m matlab/partial_information/PI_gensys.m matlab/partial_information/qzswitch.m
@@ -123,6 +123,11 @@ Copyright: 2010-2015 Alexander Meyer-Gohde
            2015-2017 Dynare Team
 License: GPL-3+
 
+Files: matlab/collapse_figures_in_tabgroup.m
+Copyright: 2023 Eduard Benet Cerda
+           2024 Dynare Team
+License: GPL-3+
+
 Files: matlab/convergence_diagnostics/raftery_lewis.m
 Copyright: 2016 Benjamin Born and Johannes Pfeifer
            2016-2017 Dynare Team
@@ -172,7 +177,7 @@ Comment: Written by Jessica Cariboni and Francesca Campolongo
 Files: matlab/+gsa/cumplot.m
        matlab/+gsa/monte_carlo_filtering.m
        matlab/+gsa/skewness.m
-       matlab/+gsa/log_trans_.m
+       matlab/+gsa/log_transform.m
        matlab/+gsa/map_calibration.m
        matlab/+gsa/map_identification.m
        matlab/+gsa/monte_carlo_filtering_analysis.m
@@ -247,7 +252,7 @@ License: BSD-2-clause
 
 Files: examples/fs2000_data.m
 Copyright: 2000-2022 Frank Schorfheide
-Copyright: 2023 Dynare Team
+           2023 Dynare Team
 License: CC-BY-SA-4.0
 
 Files: doc/*.rst doc/*.tex doc/*.svg doc/*.pdf doc/*.bib
@@ -292,28 +297,6 @@ Files: preprocessor/doc/preprocessor/*
 Copyright: 2007-2019 Dynare Team
 License: CC-BY-SA-4.0
 
-Files: contrib/jsonlab/*
-Copyright: 2011-2020 Qianqian Fang <q.fang at neu.edu>
-           2016 Bastian Bechtold
-License: GPL-3+ or BSD-3-clause
-
-Files: contrib/jsonlab/base64decode.m
-       contrib/jsonlab/base64encode.m
-       contrib/jsonlab/gzipdecode.m
-       contrib/jsonlab/gzipencode.m
-       contrib/jsonlab/zlibdecode.m
-       contrib/jsonlab/zlibencode.m
-Copyright: 2012 Kota Yamaguchi
-           2011-2020 Qianqian Fang <q.fang at neu.edu>
-License: GPL-3+ or BSD-2-clause
-
-Files: contrib/jsonlab/loadjson.m
-Copyright: 2011-2020 Qianqian Fang
-           2009 Nedialko Krouchev
-           2009 François Glineur
-           2008 Joel Feenstra
-License: GPL-3+ or BSD-2-clause or BSD-3-clause
-
 Files: contrib/ms-sbvar/utilities_dw/*
 Copyright: 1996-2011 Daniel Waggoner
 License: GPL-3+
@@ -420,32 +403,6 @@ License: BSD-2-clause
  ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  POSSIBILITY OF SUCH DAMAGE.
 
-License: BSD-3-clause
- Redistribution and use in source and binary forms, with or without modification,
- are permitted provided that the following conditions are met:
- .
- * Redistributions of source code must retain the above copyright notice, this
-   list of conditions and the following disclaimer.
- .
- * Redistributions in binary form must reproduce the above copyright notice,
-   this list of conditions and the following disclaimer in the documentation
-   and/or other materials provided with the distribution.
- .
- * Neither the name of the copyright holder nor the names of its contributors
-   may be used to endorse or promote products derived from this software without
-   specific prior written permission.
- .
- THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
- AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
- DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE
- FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
- SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
- CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
- OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
- OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
-
 License: GFDL-NIV-1.3+
  Permission is granted to copy, distribute and/or modify this document
  under the terms of the GNU Free Documentation License, Version 1.3 or

matlab/+bvar/forecast.m

@@ -119,7 +119,7 @@ OutputDirectoryName = CheckPath('graphs',M_.dname);
 dyn_graph=bvar.graph_init(sprintf('BVAR forecasts (nlags = %d)', nlags), ny, {'b-' 'g-' 'g-' 'r-' 'r-'});
 
 for i = 1:ny
-    dyn_graph=dynare_graph(dyn_graph,[ sims_no_shock_median(:, i) ...
+    dyn_graph=plot_graph(dyn_graph,[ sims_no_shock_median(:, i) ...
                         sims_no_shock_up_conf(:, i) sims_no_shock_down_conf(:, i) ...
                         sims_with_shocks_up_conf(:, i) sims_with_shocks_down_conf(:, i) ], ...
                            options_.varobs{i});
@@ -183,3 +183,31 @@ for i = 1:length(options_.varobs)
         oo_.bvar.forecast.rmse.(name) = rmse(i);
     end
 end
+
+
+function dyn_graph=plot_graph(dyn_graph,y,tit,x)
+% function plot_graph(dyn_graph, y,tit,x)
+
+if nargin < 4
+    x = (1:size(y,1))';
+end
+nplot = dyn_graph.plot_nbr + 1;
+if nplot > dyn_graph.max_nplot
+    figure('Name',dyn_graph.figure_name);
+    nplot = 1;
+end
+dyn_graph.plot_nbr = nplot;
+subplot(dyn_graph.nr,dyn_graph.nc,nplot);
+
+line_types = dyn_graph.line_types;
+line_type = line_types{1};
+for i=1:size(y,2)
+    if length(line_types) > 1
+        line_type = line_types{i};
+    end
+
+    plot(x,y(:,i),line_type);
+    hold on
+end
+title(tit);
+hold off

matlab/+gsa/map_identification.m

@@ -272,80 +272,8 @@ elseif opt_gsa.morris==3
     return
 elseif opt_gsa.morris==2   % ISKREV stuff
     return
-else  % main effects analysis
+else  
-    if itrans==0
+    error('gsa/map_identification: unsupported option morris=%u',opt_gsa.morris)
-        fsuffix = '';
-    elseif itrans==1
-        fsuffix = '_log';
-    else
-        fsuffix = '_rank';
-    end
-
-    imap=1:npT;
-
-    if isempty(lpmat0)
-        x0=lpmat(istable,:);
-    else
-        x0=[lpmat0(istable,:), lpmat(istable,:)];
-    end
-    nrun=length(istable);
-    nest=min(250,nrun);
-
-    if opt_gsa.load_ident_files==0
-        try
-            EET=load([OutputDirectoryName,'/SCREEN/',fname_,'_morris_IDE'],'SAcc','ir_cc','ic_cc');
-        catch
-            EET=[];
-        end
-        ccac = gsa.standardize_columns([mss cc ac]);
-        [pcc, dd] = eig(cov(ccac(istable,:)));
-        [latent, isort] = sort(-diag(dd));
-        latent = -latent;
-        figure, bar(latent)
-        title('Eigenvalues in PCA')
-        pcc=pcc(:,isort);
-        ccac = ccac*pcc;
-        npca = max(find(cumsum(latent)./length(latent)<0.99))+1;
-        siPCA = (EET.SAcc'*abs(pcc'))';
-        siPCA = siPCA./(max(siPCA,[],2)*ones(1,npT));
-        SAcc=zeros(size(ccac,2),npT);
-        gsa_=NaN(npca);
-        for j=1:npca %size(ccac,2),
-            if itrans==0
-                y0 = ccac(istable,j);
-            elseif itrans==1
-                y0 = gsa.log_transform(ccac(istable,j));
-            else
-                y0 = trank(ccac(istable,j));
-            end
-            if ~isempty(EET)
-                imap=find(siPCA(j,:) >= (0.1.*max(siPCA(j,:))) );
-            end
-            gsa_(j) = gsa_sdp(y0(1:nest), x0(1:nest,imap), ...
-                              2, [],[],[],0,[OutputDirectoryName,'/map_cc',fsuffix,int2str(j)], pnames);
-            SAcc(j,imap)=gsa_(j).si;
-            imap_cc{j}=imap;
-        end
-        save([OutputDirectoryName,'/map_cc',fsuffix,'.mat'],'gsa_')
-        save([OutputDirectoryName,'/',fname_,'_main_eff.mat'],'imap_cc','SAcc','ccac','-append')
-    else
-        load([OutputDirectoryName,'/',fname_,'_main_eff'],'SAcc')
-    end
-
-    hh_fig=dyn_figure(options_.nodisplay,'Name',['Identifiability indices in the ',fsuffix,' moments.']);
-    bar(sum(SAcc))
-    set(gca,'xticklabel',' ','fontsize',10,'xtick',1:npT)
-    set(gca,'xlim',[0.5 npT+0.5])
-    ydum = get(gca,'ylim');
-    set(gca,'ylim',[0 ydum(2)])
-    set(gca,'position',[0.13 0.2 0.775 0.7])
-    for ip=1:npT
-        text(ip,-0.02*(ydum(2)),bayestopt_.name{ip},'rotation',90,'HorizontalAlignment','right','interpreter','none')
-    end
-    xlabel(' ')
-    title(['Identifiability indices in the ',fsuffix,' moments.'],'interpreter','none')
-    dyn_saveas(hh_fig,[OutputDirectoryName,'/',fname_,'_ident_ALL',fsuffix],options_.nodisplay,options_.graph_format);
-    create_TeX_loader(options_,[OutputDirectoryName,'/',fname_,'_ident_ALL',fsuffix],1,['Identifiability indices in the ',fsuffix,' moments.'],['ident_ALL',fsuffix]',1)
 end
 
 function []=create_TeX_loader(options_,figpath,ifig_number,caption,label_name,scale_factor)

matlab/+gsa/run.m

@@ -13,7 +13,7 @@ function x0=run(M_,oo_,options_,bayestopt_,estim_params_,options_gsa)
 % M. Ratto (2008), Analysing DSGE Models with Global Sensitivity Analysis, 
 % Computational Economics (2008), 31, pp. 115–139
 
-% Copyright © 2008-2023 Dynare Team
+% Copyright © 2008-2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -101,6 +101,9 @@ if ~isempty(options_gsa.datafile) || isempty(bayestopt_) || options_gsa.rmse
         disp('must be specified for RMSE analysis!');
         error('Sensitivity anaysis error!')
     end
+    if isfield(options_gsa,'nobs')
+        options_.nobs=options_gsa.nobs;
+    end
     if ~isempty(options_.nobs) && length(options_.nobs)~=1
         error('dynare_sensitivity does not support recursive estimation. Please specify nobs as a scalar, not a vector.')
     end
@@ -108,9 +111,6 @@ if ~isempty(options_gsa.datafile) || isempty(bayestopt_) || options_gsa.rmse
     if isfield(options_gsa,'first_obs')
         options_.first_obs=options_gsa.first_obs;
     end
-    if isfield(options_gsa,'nobs')
-        options_.nobs=options_gsa.nobs;
-    end
     if isfield(options_gsa,'presample')
         options_.presample=options_gsa.presample;
     end
@@ -150,6 +150,11 @@ end
 
 [~,~,~,~,oo_.dr,M_.params] = dynare_resolve(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
 
+if isfield(oo_.dr,'eigval') && any(abs(oo_.dr.eigval-1)<abs(1-options_.qz_criterium)) && options_.qz_criterium<1
+    fprintf('\ngsa: The model features a unit root, but qz_criterium<1. Check whether that is intended.')
+    fprintf('\ngsa: If not, use the diffuse_filter option.\n')
+end
+
 options_gsa = set_default_option(options_gsa,'identification',0);
 if options_gsa.identification
     options_gsa.redform=0;
@@ -522,4 +527,4 @@ if options_gsa.glue
             save([OutputDirectoryName,'/',fname_,'_glue_mc.mat'], 'Out', 'Sam', 'Lik', 'Obs', 'Rem','Info', 'Exo')
         end
     end
-end
+end

matlab/+identification/run.m

@@ -736,7 +736,7 @@ if iload <=0
                 end
                 run_index = 0; % reset index
             end
-            if SampleSize > 1
+            if SampleSize > 1 && mod(iteration,3)
                 dyn_waitbar(iteration/SampleSize, h, ['MC identification checks ', int2str(iteration), '/', int2str(SampleSize)]);
             end
         end

matlab/+identification/simulated_moment_uncertainty.m

@@ -12,7 +12,7 @@ function [cmm, mm] = simulated_moment_uncertainty(indx, periods, replic,options_
 %   - cmm:      [n_moments by n_moments] covariance matrix of simulated moments
 %   - mm:       [n_moments by replic] matrix of moments
 
-% Copyright © 2009-2018 Dynare Team
+% Copyright © 2009-2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -60,15 +60,10 @@ end
 
 oo_.dr=set_state_space(oo_.dr,M_);
 
-
 if options_.logged_steady_state %if steady state was previously logged, undo this
     oo_.dr.ys=exp(oo_.dr.ys);
     oo_.steady_state=exp(oo_.steady_state);
     options_.logged_steady_state=0;
-    logged_steady_state_indicator=1;
-    evalin('base','options_.logged_steady_state=0;')
-else
-    logged_steady_state_indicator=0;
 end
 
 [oo_.dr,info,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
@@ -92,7 +87,6 @@ else
     end
 end
 
-
 for j=1:replic
     [ys, oo_.exo_simul] = simult(y0,oo_.dr,M_,options_);%do simulation
     oo_=disp_moments(ys, options_.varobs,M_,options_,oo_); %get moments
@@ -106,8 +100,5 @@ for j=1:replic
 end
 dyn_waitbar_close(h);
 
-if logged_steady_state_indicator
-    evalin('base','options_.logged_steady_state=1;') %reset base workspace option to conform to base oo_
-end
 cmm = cov(mm');
 disp('Simulated moment uncertainty ... done!')

matlab/+occbin/DSGE_smoother.m

@@ -479,7 +479,7 @@ if (~is_changed || occbin_smoother_debug) && nargin==12
                         fprintf(fidTeX,'\\begin{figure}[H]\n');
                         fprintf(fidTeX,'\\centering \n');
                         fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_smoothedshocks_occbin%s}\n',options_.figures.textwidth*min(j1/3,1),[GraphDirectoryName '/' M_.fname],int2str(ifig)); % don't use filesep as it will create issues with LaTeX on Windows
-                        fprintf(fidTeX,'\\caption{Check plots.}');
+                        fprintf(fidTeX,'\\caption{OccBin smoothed shocks.}');
                         fprintf(fidTeX,'\\label{Fig:smoothedshocks_occbin:%s}\n',int2str(ifig));
                         fprintf(fidTeX,'\\end{figure}\n');
                         fprintf(fidTeX,' \n');
@@ -488,7 +488,7 @@ if (~is_changed || occbin_smoother_debug) && nargin==12
             end
         end
 
-            if mod(j1,9)~=0 && j==M_.exo_nbr
+        if mod(j1,9)~=0 && j==M_.exo_nbr
                 annotation('textbox', [0.1,0,0.35,0.05],'String', 'Linear','Color','Blue','horizontalalignment','center','interpreter','none');
                 annotation('textbox', [0.55,0,0.35,0.05],'String', 'Piecewise','Color','Red','horizontalalignment','center','interpreter','none');
                 dyn_saveas(hh_fig,[GraphDirectoryName filesep M_.fname,'_smoothedshocks_occbin',int2str(ifig)],options_.nodisplay,options_.graph_format);
@@ -497,7 +497,7 @@ if (~is_changed || occbin_smoother_debug) && nargin==12
                 fprintf(fidTeX,'\\begin{figure}[H]\n');
                 fprintf(fidTeX,'\\centering \n');
                 fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_smoothedshocks_occbin%s}\n',options_.figures.textwidth*min(j1/3,1),[GraphDirectoryName '/' M_.fname],int2str(ifig)); % don't use filesep as it will create issues with LaTeX on Windows
-                fprintf(fidTeX,'\\caption{Check plots.}');
+                fprintf(fidTeX,'\\caption{OccBin smoothed shocks.}');
                 fprintf(fidTeX,'\\label{Fig:smoothedshocks_occbin:%s}\n',int2str(ifig));
                 fprintf(fidTeX,'\\end{figure}\n');
                 fprintf(fidTeX,' \n');
@@ -505,6 +505,6 @@ if (~is_changed || occbin_smoother_debug) && nargin==12
         end
         if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
             fclose(fidTeX);
+        end
     end
 end
-end

matlab/+occbin/forecast.m

@@ -1,143 +1,169 @@
-function [oo_, error_flag] = forecast(options_,M_,oo_,forecast) %,hist_period)
+function [forecast, error_flag] = forecast(options_,M_,dr,endo_steady_state,exo_steady_state,exo_det_steady_state,forecast_horizon)
-%function oo_ = forecast(options_,M_,oo_,forecast)
+% [forecast, error_flag] = forecast(options_,M_,dr,endo_steady_state,exo_steady_state,exo_det_steady_state,forecast_horizon)
-% forecast
+% Occbin forecasts
+%
+% INPUTS
+% - options_                [structure]     Matlab's structure describing the current options
+% - M_                      [structure]     Matlab's structure describing the model
+% - dr_in                   [structure]     model information structure
+% - endo_steady_state       [double]        steady state value for endogenous variables
+% - exo_steady_state        [double]        steady state value for exogenous variables
+% - exo_det_steady_state    [double]        steady state value for exogenous deterministic variables                                    
+% - forecast_horizon        [integer]       forecast horizon
+%
+% OUTPUTS
+% - forecast                [structure]     forecast results
+% - error_flag              [integer]       error code
+%
+% SPECIAL REQUIREMENTS
+%   none.
+
+% Copyright © 2022-2023 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 opts = options_.occbin.forecast;
 
 options_.occbin.simul.maxit = opts.maxit;
 options_.occbin.simul.check_ahead_periods = opts.check_ahead_periods;
-options_.occbin.simul.periods = forecast;
+options_.occbin.simul.periods = forecast_horizon;
-SHOCKS0 = opts.SHOCKS0;
+shocks_input = opts.SHOCKS0;
-if ~isempty(SHOCKS0)
+
-    if iscell(SHOCKS0)
+if ~isempty(shocks_input)
-        for j=1:length(SHOCKS0)
+    n_shocks=size(shocks_input,1);
-            sname = SHOCKS0{j}{1};
+    if iscell(shocks_input)
-            inds(j)=strmatch(sname,M_.exo_names);
+        inds=NaN(n_shocks,1);
-            SHOCKS1(j,:)=SHOCKS0{j}{2};
+        periods=length(shocks_input{1}{2});
+        shock_mat=NaN(n_shocks,periods);
+        for j=1:n_shocks
+            exo_pos=strmatch(shocks_input{j}{1},M_.exo_names,'exact');
+            if isempty(exo_pos)
+        	    error('occbin.forecast: unknown exogenous shock %s',shocks_input{j}{1})
+            else
+                inds(j)=exo_pos;
+            end
+            if length(shocks_input{j}{2})~=periods
+        	    error('occbin.forecast: unknown exogenous shock %s',shocks_input{j}{1})
+            else
+                shock_mat(j,:)=shocks_input{j}{2};
+            end
         end
-    elseif isreal(SHOCKS0)
+    elseif isreal(shocks_input)
-        SHOCKS1=SHOCKS0;
+        shock_mat=shocks_input;
-        inds = 1:M_.exo_nbr;
+        inds = (1:M_.exo_nbr)';
     end
 end
 
+options_.occbin.simul.endo_init = M_.endo_histval(:,1)-endo_steady_state; %initial condition
+options_.occbin.simul.init_regime = opts.frcst_regimes;
+options_.occbin.simul.init_binding_indicator = [];
+
+shocks_base = zeros(forecast_horizon,M_.exo_nbr);
+if ~isempty(shocks_input)
+    for j=1:n_shocks
+        shocks_base(:,inds(j))=shock_mat(j,:);
+    end
+end
 
 if opts.replic
     h = dyn_waitbar(0,'Please wait occbin forecast replic ...');
     ishock = find(sqrt(diag((M_.Sigma_e))));
+    options_.occbin.simul.exo_pos=ishock;
     effective_exo_nbr=  length(ishock);
-    effective_exo_names = M_.exo_names(ishock);
+    effective_Sigma_e = M_.Sigma_e(ishock,ishock);  % does not take heteroskedastic shocks into account
-    effective_Sigma_e = M_.Sigma_e(ishock,ishock);
     [U,S] = svd(effective_Sigma_e);
+    % draw random shocks
     if opts.qmc
         opts.replic =2^(round(log2(opts.replic+1)))-1;
-        SHOCKS_ant =   qmc_sequence(forecast*effective_exo_nbr, int64(1), 1, opts.replic)';
+        SHOCKS_add = qmc_sequence(forecast_horizon*effective_exo_nbr, int64(1), 1, opts.replic);
     else
-        SHOCKS_ant = randn(forecast*effective_exo_nbr,opts.replic)';
+        SHOCKS_add = randn(forecast_horizon*effective_exo_nbr,opts.replic);
     end
-    zlin0=zeros(forecast,M_.endo_nbr,opts.replic);
+    SHOCKS_add=reshape(SHOCKS_add,effective_exo_nbr,forecast_horizon,opts.replic);
-    zpiece0=zeros(forecast,M_.endo_nbr,opts.replic);
+    z.linear=NaN(forecast_horizon,M_.endo_nbr,opts.replic);
+    z.piecewise=NaN(forecast_horizon,M_.endo_nbr,opts.replic);
+    error_flag=true(opts.replic,1);
+    simul_SHOCKS=NaN(forecast_horizon,M_.exo_nbr,opts.replic);
     for iter=1:opts.replic
-
+        options_.occbin.simul.SHOCKS = shocks_base+transpose(U*sqrt(S)*SHOCKS_add(:,:,iter));
-        if ~isempty(SHOCKS0)
-            for j=1:length(SHOCKS0)
-                SHOCKS(:,inds(j))=SHOCKS1(j,:);
-            end
-        end
-
-        error_flagx=1;
-        % while error_flagx,
-        % SHOCKS=transpose(sqrt(diag(diag(effective_Sigma_e)))*(reshape(SHOCKS_ant(iter,:),forecast,effective_exo_nbr))');
-        SHOCKS=transpose(U*sqrt(S)*(reshape(SHOCKS_ant(iter,:),forecast,effective_exo_nbr))');
-        %             SHOCKS=transpose(U*sqrt(S)*randn(forecast,M_.exo_nbr)');   %realized shocks
-        options_.occbin.simul.endo_init = M_.endo_histval(:,1)-oo_.steady_state;
-        options_.occbin.simul.init_regime = opts.frcst_regimes;
-        options_.occbin.simul.init_binding_indicator = [];
-        options_.occbin.simul.exo_pos=ishock;
-        options_.occbin.simul.SHOCKS = SHOCKS;
         options_.occbin.simul.waitbar=0;
-        [~, out] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
+        [~, out] = occbin.solver(M_,options_,dr,endo_steady_state,exo_steady_state,exo_det_steady_state);
-        zlin0(:,:,iter)=out.linear;
-        zpiece0(:,:,iter)=out.piecewise;
-        ys=out.ys;
-        frcst_regime_history(iter,:)=out.regime_history;
         error_flag(iter)=out.error_flag;
-        error_flagx = error_flag(iter);
+        if ~error_flag(iter)
-        % end
+            z.linear(:,:,iter)=out.linear;
-        simul_SHOCKS(:,:,iter) = SHOCKS;
+            z.piecewise(:,:,iter)=out.piecewise;
-
+            frcst_regime_history(iter,:)=out.regime_history;
-        if error_flag(iter) && debug_flag
+            error_flag(iter)=out.error_flag;
-            %             display('no solution')
+            simul_SHOCKS(:,:,iter) = shocks_base;
-
+        else
-%             keyboard;
+            if options_.debug
-            save no_solution SHOCKS zlin0 zpiece0 iter frcst_regime_history
+                save('Occbin_forecast_debug','simul_SHOCKS','z','iter','frcst_regime_history','error_flag','out','shocks_base')
+            end
         end
         dyn_waitbar(iter/opts.replic,h,['OccBin MC forecast replic ',int2str(iter),'/',int2str(opts.replic)])
     end
     dyn_waitbar_close(h);
-    save temp zlin0 zpiece0 simul_SHOCKS error_flag
+    if options_.debug
+         save('Occbin_forecast_debug','simul_SHOCKS','z','iter','frcst_regime_history','error_flag')
+    end
     inx=find(error_flag==0);
-    zlin0=zlin0(:,:,inx);
+    z.linear=z.linear(:,:,inx);
-    zpiece0=zpiece0(:,:,inx);
+    z.piecewise=z.piecewise(:,:,inx);
-    zlin   = mean(zlin0,3);
+    z.min.piecewise = min(z.piecewise,[],3);
-    zpiece = mean(zpiece0,3);
+    z.max.piecewise = max(z.piecewise,[],3);
-    zpiecemin = min(zpiece0,[],3);
+    z.min.linear = min(z.linear,[],3);
-    zpiecemax = max(zpiece0,[],3);
+    z.max.linear = max(z.linear,[],3);
-    zlinmin = min(zlin0,[],3);
+    
-    zlinmax = max(zlin0,[],3);
+    field_names={'linear','piecewise'};
-
+    post_mean=NaN(forecast_horizon,1);
-    for i=1:M_.endo_nbr
+    post_median=NaN(forecast_horizon,1);
-        for j=1:forecast
+    post_var=NaN(forecast_horizon,1);
-            [post_mean(j,1), post_median(j,1), post_var(j,1), hpd_interval(j,:), post_deciles(j,:)] = posterior_moments(squeeze(zlin0(j,i,:)),options_.forecasts.conf_sig);
+    hpd_interval=NaN(forecast_horizon,2);
+    post_deciles=NaN(forecast_horizon,9);
+    for field_iter=1:2
+        for i=1:M_.endo_nbr
+            for j=1:forecast_horizon
+                [post_mean(j,1), post_median(j,1), post_var(j,1), hpd_interval(j,:), post_deciles(j,:)] = posterior_moments(squeeze(z.(field_names{field_iter})(j,i,:)),options_.forecasts.conf_sig);
+            end
+            forecast.(field_names{field_iter}).Mean.(M_.endo_names{i})=post_mean;
+            forecast.(field_names{field_iter}).Median.(M_.endo_names{i})=post_median;
+            forecast.(field_names{field_iter}).Var.(M_.endo_names{i})=post_var;
+            forecast.(field_names{field_iter}).HPDinf.(M_.endo_names{i})=hpd_interval(:,1);
+            forecast.(field_names{field_iter}).HPDsup.(M_.endo_names{i})=hpd_interval(:,2);
+            forecast.(field_names{field_iter}).Deciles.(M_.endo_names{i})=post_deciles;
+            forecast.(field_names{field_iter}).Min.(M_.endo_names{i})=z.min.(field_names{field_iter})(:,i);
+            forecast.(field_names{field_iter}).Max.(M_.endo_names{i})=z.max.(field_names{field_iter})(:,i);
         end
-        oo_.occbin.linear_forecast.Mean.(M_.endo_names{i})=post_mean;
-        oo_.occbin.linear_forecast.Median.(M_.endo_names{i})=post_median;
-        oo_.occbin.linear_forecast.Var.(M_.endo_names{i})=post_var;
-        oo_.occbin.linear_forecast.HPDinf.(M_.endo_names{i})=hpd_interval(:,1);
-        oo_.occbin.linear_forecast.HPDsup.(M_.endo_names{i})=hpd_interval(:,2);
-        oo_.occbin.linear_forecast.Deciles.(M_.endo_names{i})=post_deciles;
-        oo_.occbin.linear_forecast.Min.(M_.endo_names{i})=zlinmin(:,i);
-        oo_.occbin.linear_forecast.Max.(M_.endo_names{i})=zlinmax(:,i);
-        for j=1:forecast
-            [post_mean(j,1), post_median(j,1), post_var(j,1), hpd_interval(j,:), post_deciles(j,:)] = posterior_moments(squeeze(zpiece0(j,i,:)),options_.forecasts.conf_sig);
-        end
-        oo_.occbin.forecast.Mean.(M_.endo_names{i})=post_mean;
-        oo_.occbin.forecast.Median.(M_.endo_names{i})=post_median;
-        oo_.occbin.forecast.Var.(M_.endo_names{i})=post_var;
-        oo_.occbin.forecast.HPDinf.(M_.endo_names{i})=hpd_interval(:,1);
-        oo_.occbin.forecast.HPDsup.(M_.endo_names{i})=hpd_interval(:,2);
-        oo_.occbin.forecast.Deciles.(M_.endo_names{i})=post_deciles;
-        oo_.occbin.forecast.Min.(M_.endo_names{i})=zpiecemin(:,i);
-        oo_.occbin.forecast.Max.(M_.endo_names{i})=zpiecemax(:,i);
-        %   eval([M_.endo_names{i},'_ss=zdatass(i);']);
     end
-
 else
-    SHOCKS = zeros(forecast,M_.exo_nbr);
-    if ~isempty(SHOCKS0)
-        for j=1:length(SHOCKS0)
-            SHOCKS(:,inds(j))=SHOCKS1(j,:);
-        end
-    end
-    options_.occbin.simul.endo_init = M_.endo_histval(:,1)-oo_.steady_state;
-    options_.occbin.simul.init_regime = opts.frcst_regimes;
-    options_.occbin.simul.init_violvecbool = [];
     options_.occbin.simul.irfshock = M_.exo_names;
-    options_.occbin.simul.SHOCKS = SHOCKS;
+    options_.occbin.simul.SHOCKS = shocks_base;
-    [~, out] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
+    [~, out] = occbin.solver(M_,options_,dr,endo_steady_state,exo_steady_state,exo_det_steady_state);
+    error_flag=out.error_flag;
+    if out.error_flag
+        fprintf('occbin.forecast: forecast simulation failed.')
+        return;
+    end
 
-    zlin=out.linear;
-    zpiece=out.piecewise;
     frcst_regime_history=out.regime_history;
     error_flag=out.error_flag;
     for i=1:M_.endo_nbr
-        oo_.occbin.linear_forecast.Mean.(M_.endo_names{i})= zlin(:,i);
+        forecast.linear.Mean.(M_.endo_names{i})= out.linear(:,i);
-        oo_.occbin.forecast.Mean.(M_.endo_names{i})= zpiece(:,i);
+        forecast.piecewise.Mean.(M_.endo_names{i})= out.piecewise(:,i);
-        oo_.occbin.forecast.HPDinf.(M_.endo_names{i})= nan;
-        oo_.occbin.forecast.HPDsup.(M_.endo_names{i})= nan;
     end
-
 end
 
-oo_.occbin.forecast.regimes=frcst_regime_history;
+forecast.regimes=frcst_regime_history;
-

matlab/+occbin/irf.m

@@ -1,140 +1,124 @@
-function [oo_] = irf(M_,oo_,options_)
+function irfs = irf(M_,oo_,options_)
+% irfs = irf(M_,oo_,options_)
+% Calls a minimizer
+%
+% INPUTS
+% - M_                  [structure]     Matlab's structure describing the model
+% - oo_                 [structure]     Matlab's structure containing the results
+% - options_            [structure]     Matlab's structure describing the current options
+%
+% OUTPUTS
+% - irfs                [structure]     IRF results
+%
+% SPECIAL REQUIREMENTS
+%   none.
+%
+%
+% Copyright © 2022-2023 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 shocknames = options_.occbin.irf.exo_names;
-shocksigns = options_.occbin.irf.shocksigns;
+shocksigns = options_.occbin.irf.shocksigns; %'pos','neg'
 shocksize = options_.occbin.irf.shocksize;
-t0 = options_.occbin.irf.t0;
+t_0 = options_.occbin.irf.t0;
+
+%% set simulation options based on IRF options
 options_.occbin.simul.init_regime = options_.occbin.irf.init_regime;
 options_.occbin.simul.check_ahead_periods = options_.occbin.irf.check_ahead_periods;
 options_.occbin.simul.maxit = options_.occbin.irf.maxit;
 options_.occbin.simul.periods = options_.irf;
 
-% Run inital conditions + other shocks
+%% Run initial conditions + other shocks
-if t0 == 0
+if t_0 == 0
-    shocks0= zeros(options_.occbin.simul.periods+1,M_.exo_nbr);
+    shocks_base = zeros(options_.occbin.simul.periods+1,M_.exo_nbr);
     options_.occbin.simul.endo_init = [];
 else
-    % girf conditional to smoothed states in t0 and shocks in t0+1
+    if ~isfield(oo_.occbin,'smoother')
-    shocks0= [oo_.occbin.smoother.etahat(:,t0+1)'; zeros(options_.occbin.simul.periods,M_.exo_nbr)];
+        error('occbin.irfs: smoother must be run before requesting GIRFs based on smoothed results')
-    options_.occbin.simul.SHOCKS=shocks0;
+    end
-    options_.occbin.simul.endo_init = oo_.occbin.smoother.alphahat(oo_.dr.inv_order_var,t0);
+    % GIRF conditional on smoothed states in t_0 and shocks in t_0+1
+    shocks_base= [oo_.occbin.smoother.etahat(:,t_0+1)'; zeros(options_.occbin.simul.periods,M_.exo_nbr)];
+    options_.occbin.simul.SHOCKS=shocks_base;
+    options_.occbin.simul.endo_init = oo_.occbin.smoother.alphahat(oo_.dr.inv_order_var,t_0);
 end
-options_.occbin.simul.SHOCKS=shocks0;
+options_.occbin.simul.SHOCKS=shocks_base;
-[~, out0] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
-zlin0 = out0.linear;
-zpiece0 = out0.piecewise;
 
-% Select shocks of interest
+[~, out_base] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
-jexo_all =zeros(size(shocknames,1),1);
+if out_base.error_flag
+    error('occbin.irfs: could not compute the solution')
+end
 
+irfs.linear = struct();
+irfs.piecewise = struct();
+
+% Get indices of shocks of interest
+exo_index =zeros(size(shocknames,1),1);
 for i=1:length(shocknames)
-    jexo_all(i) = strmatch(shocknames{i},M_.exo_names,'exact');
+    exo_index(i) = strmatch(shocknames{i},M_.exo_names,'exact');
 end
 
-oo_.occbin.linear_irfs = struct();
+% cs=get_lower_cholesky_covariance(M_.Sigma_e,options_.add_tiny_number_to_cholesky);
-oo_.occbin.irfs = struct();
+% irf_shocks_indx = getIrfShocksIndx(M_, options_);
+
 % Set shock size
 if isempty(shocksize)
-    %         if isfield(oo_.posterior_mode.shocks_std,M_.exo_names{jexo})
+    shocksize = sqrt(diag(M_.Sigma_e(exo_index,exo_index)));
-    shocksize = sqrt(diag(M_.Sigma_e(jexo_all,jexo_all))); %oo_.posterior_mode.shocks_std.(M_.exo_names{jexo});
     if any(shocksize < 1.e-9)
         shocksize(shocksize < 1.e-9) = 0.01;
     end
 end
+
 if numel(shocksize)==1
     shocksize=repmat(shocksize,[length(shocknames),1]);
 end
 
 % Run IRFs
-for counter = 1:length(jexo_all)
+for sign_iter=1:length(shocksigns)
-  
+    for IRF_counter = 1:length(exo_index)
-    jexo = jexo_all(counter);
+        jexo = exo_index(IRF_counter);
-    
+        if ~options_.noprint && options_.debug
-    if ~options_.noprint
+            fprintf('occbin.irf: Producing GIRFs for shock %s. Simulation %d out of %d. \n',M_.exo_names{jexo},IRF_counter,size(exo_index,1));
-        fprintf('Producing GIRFs for shock %s. Simulation %d out of %d. \n',M_.exo_names{jexo},counter,size(jexo_all,1));
+        end
-    end 
+        shocks1=shocks_base;
-    
+        if ismember('pos',shocksigns{sign_iter})
-    if ismember('pos',shocksigns)
+            shocks1(1,jexo)=shocks_base(1,jexo)+shocksize(IRF_counter);
-        % (+) shock
+        elseif ismember('neg',shocksigns{sign_iter})
-        shocks1=shocks0;
+            shocks1(1,jexo)=shocks_base(1,jexo)-shocksize(IRF_counter);
-        shocks1(1,jexo)=shocks0(1,jexo)+shocksize(counter);
+        end
-        if t0 == 0
+        options_.occbin.simul.SHOCKS=shocks1;
-            options_.occbin.simul.SHOCKS=shocks1;
+        if t_0 == 0
             options_.occbin.simul.endo_init = [];
-            [~, out_pos] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
         else
-            options_.occbin.simul.SHOCKS=shocks1;
+            options_.occbin.simul.endo_init = oo_.occbin.smoother.alphahat(oo_.dr.inv_order_var,t_0);
-            options_.occbin.simul.endo_init = oo_.occbin.smoother.alphahat(oo_.dr.inv_order_var,t0);
-            [~, out_pos] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
         end
-        if out_pos.error_flag
+        [~, out_sim] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
-            warning('Occbin error.')
+        if out_sim.error_flag
-            return
+            warning('occbin.irfs: simulation failed')
+            skip
         end
-        zlin_pos = out_pos.linear;
-        zpiece_pos = out_pos.piecewise;
         % Substract inital conditions + other shocks
-        zlin_pos_diff       = zlin_pos-zlin0;
+        zdiff.linear.(shocksigns{sign_iter}) = out_sim.linear-out_base.linear;
-        zpiece_pos_diff      = zpiece_pos-zpiece0;
+        zdiff.piecewise.(shocksigns{sign_iter}) = out_sim.piecewise-out_base.piecewise;
-    end
+
-    
+        for j_endo=1:M_.endo_nbr
-    if ismember('neg',shocksigns)
+            if ismember('pos',shocksigns)
-        % (-) shock
+                irfs.piecewise.([M_.endo_names{j_endo} '_' M_.exo_names{jexo} '_' shocksigns{sign_iter}])   = zdiff.piecewise.(shocksigns{sign_iter})(:,j_endo);
-        shocks_1=shocks0;
+                irfs.linear.([M_.endo_names{j_endo} '_' M_.exo_names{jexo} '_' shocksigns{sign_iter}])   = zdiff.linear.(shocksigns{sign_iter})(:,j_endo);
-        shocks_1(1,jexo)=shocks0(1,jexo)-shocksize(counter);
+            end
-        if t0 == 0
-            options_.occbin.simul.SHOCKS=shocks_1;
-            options_.occbin.simul.endo_init = [];
-            [~, out_neg] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
-        else
-            options_.occbin.simul.SHOCKS=shocks_1;
-            options_.occbin.simul.endo_init = oo_.occbin.smoother.alphahat(oo_.dr.inv_order_var,t0);
-            [~, out_neg] = occbin.solver(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
         end
-        if out_neg.error_flag
-            warning('Occbin error.')
-            return
-        end
-        zlin_neg    = out_neg.linear;
-        zpiece_neg  = out_neg.piecewise;
-        zlin_neg_diff    = zlin_neg-zlin0;
-        zpiece_neg_diff  = zpiece_neg-zpiece0;
     end
-
+end
-    % Save
-    if ~isfield(oo_.occbin,'linear_irfs')
-        oo_.occbin.linear_irfs = struct();
-    end
-    if ~isfield(oo_.occbin,'irfs')
-        oo_.occbin.irfs = struct();
-    end
-
-    for jendo=1:M_.endo_nbr
-%         oo_.occbin.irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '1'])   = zpiece_pos (:,jendo);
-%         oo_.occbin.irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_1'])  = zpiece_neg (:,jendo);
-%         oo_.occbin.linear_irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '1'])   = zlin_pos (:,jendo);
-%         oo_.occbin.linear_irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_1'])  = zlin_neg(:,jendo);
-   
-        if ismember('pos',shocksigns)
-        oo_.occbin.irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_pos'])   = zpiece_pos_diff (:,jendo);
-        oo_.occbin.linear_irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_pos'])   = zlin_pos_diff (:,jendo);
-        end
-        
-        if ismember('neg',shocksigns)
-            oo_.occbin.irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_neg'])  = zpiece_neg_diff (:,jendo);
-            oo_.occbin.linear_irfs.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_neg'])  = zlin_neg_diff (:,jendo);
-        end
-        
-% %         
-%         oo_.occbin.irfs0.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '1'])   = zpiece0(:,jendo);
-%         oo_.occbin.linear_irfs0.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '1'])  = zlin0(:,jendo);
-%         oo_.occbin.irfs0.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_1'])   = zpiece0(:,jendo);
-%         oo_.occbin.linear_irfs0.([M_.endo_names{jendo} '_' M_.exo_names{jexo} '_1'])  = zlin0(:,jendo);
-
-    end
-end
-
-
-
-end
-

matlab/+occbin/plot_irfs.m

@@ -1,34 +1,74 @@
-function  plot_irfs(M_,oo_,options_,irf3,irf4)
+function plot_irfs(M_,irfs,options_,var_list)
+% plot_irfs(M_,irfs,options_,var_list)
+%
+% INPUTS
+% - M_                      [structure]     Matlab's structure describing the model
+% - irfs                    [structure]     IRF results
+% - options_                [structure]     Matlab's structure describing the current options
+% - var_list                [character array]  list of endogenous variables specified
+%
+% OUTPUTS
+%   none
+%
+% SPECIAL REQUIREMENTS
+%   none.
 
-shocknames = options_.occbin.plot_irf.exo_names;
+% Copyright © 2022-2023 Dynare Team
-simulname = options_.occbin.plot_irf.simulname;
+%
-if isempty(simulname)
+% This file is part of Dynare.
-   simulname_ = simulname;
+%
-else 
+% Dynare is free software: you can redistribute it and/or modify
-   simulname_ = [ simulname '_' ];
+% 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 <https://www.gnu.org/licenses/>.
+
+if (isfield(options_,'irf_shocks')==0)
+    shock_names  = M_.exo_names;
+else
+    shock_names  = options_.irf_shocks;
 end
-vars_irf = options_.occbin.plot_irf.endo_names;
+
-endo_names_long = options_.occbin.plot_irf.endo_names_long;
+simul_name = options_.occbin.plot_irf.simulname;
+if isempty(simul_name)
+    save_name = simul_name;
+else
+    save_name = [ simul_name '_' ];
+end
+
+if isempty(var_list)
+    var_list = M_.endo_names(1:M_.orig_endo_nbr);
+end
+
+[i_var, ~, index_uniques] = varlist_indices(var_list, M_.endo_names);
+vars_irf=var_list(index_uniques);
+
 endo_scaling_factor = options_.occbin.plot_irf.endo_scaling_factor;
-length_irf = options_.occbin.plot_irf.tplot;
+length_irf = options_.irf;
-if isempty(length_irf)
-    length_irf = options_.irf;
-end
 
-irflocation_lin   = oo_.occbin.linear_irfs;
+steps_irf = 1;
-irflocation_piece = oo_.occbin.irfs;
-
-
-steps_irf = 1; 
-warning('off','all')
 
 DirectoryName = CheckPath('graphs',M_.dname);
+latexFolder = CheckPath('latex',M_.dname);
+if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
+    fidTeX = fopen([latexFolder '/' M_.fname '_occbin_irfs.tex'],'w');
+    fprintf(fidTeX,'%% TeX eps-loader file generated by occbin.plot_irfs.m (Dynare).\n');
+    fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
+    fprintf(fidTeX,' \n');
+end
 
 iexo=[];
-for i=1:size(shocknames,1)
+for var_iter=1:size(shock_names,1)
-    itemp = strmatch(shocknames{i},M_.exo_names,'exact');
+    itemp = strmatch(shock_names{var_iter},M_.exo_names,'exact');
     if isempty(itemp)
-        error(['Shock ',shocknames{i},' is not defined!'])
+        error(['Shock ',shock_names{var_iter},' is not defined!'])
     else
         iexo=[iexo, itemp];
     end
@@ -38,104 +78,102 @@ ncols     = options_.occbin.plot_irf.ncols;
 nrows     = options_.occbin.plot_irf.nrows;
 npan      = ncols*nrows;
 
-plot_grid = options_.occbin.plot_irf.grid;
+shocksigns = options_.occbin.plot_irf.shocksigns;
-shocksigns = options_.occbin.plot_irf.shocksigns; 
+threshold = options_.impulse_responses.plot_threshold;
-threshold = options_.occbin.plot_irf.threshold; 
-
-% Add steady_state
-if options_.occbin.plot_irf.add_steadystate
-    add_stst = options_.occbin.plot_irf.add_steadystate;
-else
-    add_stst = 0;
-end 
-for sss = 1:numel(shocksigns)
-    
-    shocksign = shocksigns{sss};
-    
-    for j=1:size(shocknames,1)
-        %shocknames = M_.exo_names{j};
 
+for shock_sign_iter = 1:numel(shocksigns)
+    shocksign = shocksigns{shock_sign_iter};
+    if strcmp(shocksign,'pos')
+        plot_title_sign='positive';
+    elseif strcmp(shocksign,'neg')
+        plot_title_sign='negative';
+    else
+       error('Unknown shock sign %s',shocksign);
+    end
+        
+    for shock_iter=1:size(shock_names,1)
         j1   = 0;
         isub = 0;
         ifig = 0;
-
         % Variables
-        % ----------------------
+        for var_iter = 1:length(vars_irf)
-        for i = 1:length(vars_irf)
-
             j1=j1+1;
             if mod(j1,npan)==1
-               % vector corresponds to [left bottom width height]. 680 and 678 for the left and bottom elements correspond to the default values used by MATLAB while creating a figure and width, .
+                % vector corresponds to [left bottom width height]. 680 and 678 for the left and bottom elements correspond to the default values used by MATLAB while creating a figure and width, .
-                hfig = dyn_figure(options_.nodisplay,'name',['OccbinIRFs ' shocknames{j} ' ' simulname ' ' shocksign],'PaperPositionMode', 'auto','PaperType','A4','PaperOrientation','portrait','renderermode','auto','position',[10 10 950 650]);
+                screensize = get( groot, 'Screensize' );
+                hfig = dyn_figure(options_.nodisplay,'name',['OccBin IRFs to ' plot_title_sign ' ' shock_names{shock_iter} ' shock ' simul_name],'OuterPosition' ,[50 50 min(1000,screensize(3)-50) min(750,screensize(4)-50)]);
                 ifig=ifig+1;
                 isub=0;
             end
-            isub=isub+1;      
+            isub=isub+1;
-            
             if isempty(endo_scaling_factor)
                 exofactor = 1;
-            else               
-                exofactor = endo_scaling_factor{i};
-            end 
-
-            subplot(nrows,ncols,isub)            
-            irf_field   = strcat(vars_irf{i,1},'_',shocknames{j},'_',shocksign);
-            
-            irfvalues   = irflocation_lin.(irf_field);
-            if add_stst
-                irfvalues = irfvalues + get_mean(vars_irf{i,1});
-            end          
-            
-            irfvalues(abs(irfvalues) <threshold) = 0;
-            
-            plot(irfvalues(1:steps_irf:length_irf)*exofactor,'linewidth',2);
-            hold on
-
-            irfvalues   = irflocation_piece.(irf_field);  
-            if add_stst
-                irfvalues = irfvalues + get_mean(vars_irf{i,1});
-            end
-            irfvalues(abs(irfvalues) <threshold) = 0;
-            plot(irfvalues(1:steps_irf:length_irf)*exofactor,'r--','linewidth',2);
-            
-            hold on
-            plot(irfvalues(1:steps_irf:length_irf)*0,'k-','linewidth',1.5);
-            % Optional additional IRFs
-            if nargin > 10
-                  irfvalues   = irf3.(irf_field) ;
-                  irfvalues(abs(irfvalues) <threshold) = 0;
-                  plot(irfvalues(1:steps_irf:length_irf)*exofactor,'k:','linewidth',2);
-            end
-            if nargin > 11
-                  irfvalues   = irf4.(irf_field)   ;   
-                  irfvalues(abs(irfvalues) <threshold) = 0;
-                  plot(irfvalues(1:steps_irf:length_irf)*exofactor,'g-.','linewidth',2);
-            end
-      
-
-            if plot_grid
-                grid on
-            end
-
-            xlim([1 (length_irf/steps_irf)]);
-            
-            % title
-            if isempty(endo_names_long)
-                title(regexprep(vars_irf{i},'_',' '))
             else
-                title(endo_names_long{i})
+                exofactor = endo_scaling_factor{var_iter};
             end
 
+            subplot(nrows,ncols,isub)
+            irf_field   = strcat(vars_irf{var_iter},'_',shock_names{shock_iter},'_',shocksign);
+
+            %%linear IRFs
+            if ~isfield(irfs.linear,irf_field)
+                warning('occbin.plot_irfs: no linear IRF for %s to %s detected',vars_irf{var_iter,1},shock_names{shock_iter})
+            else
+                irfvalues   = irfs.linear.(irf_field);
+                irfvalues(abs(irfvalues) <threshold) = 0;
+                if options_.occbin.plot_irf.add_steadystate
+                    irfvalues = irfvalues + get_mean(vars_irf{var_iter,1});
+                end
+                max_irf_length_1=min(length_irf,length(irfvalues));
+                plot(irfvalues(1:steps_irf:max_irf_length_1)*exofactor,'linewidth',2);
+            end
+            hold on
+
+            if ~isfield(irfs.piecewise,irf_field)
+                warning('occbin.plot_irfs: no piecewise linear IRF for %s to %s detected',vars_irf{var_iter,1},shock_names{shock_iter})
+            else
+                irfvalues   = irfs.piecewise.(irf_field);
+                irfvalues(abs(irfvalues) <threshold) = 0;
+                if options_.occbin.plot_irf.add_steadystate
+                    irfvalues = irfvalues + get_mean(vars_irf{var_iter,1});
+                end
+                max_irf_length_2=min(length_irf,length(irfvalues));
+                plot(irfvalues(1:steps_irf:max_irf_length_2)*exofactor,'r--','linewidth',2);
+            end
+
+            plot(irfvalues(1:steps_irf:max(max_irf_length_1,max_irf_length_2))*0,'k-','linewidth',1.5);
+
+            if options_.occbin.plot_irf.grid
+                grid on
+            end     
+            xlim([1 max(max_irf_length_1,max_irf_length_2)]);
+            if options_.TeX
+                title(['$' M_.endo_names_tex{i_var(var_iter)}, '$'],'Interpreter','latex')
+            else
+                title(M_.endo_names{i_var(var_iter)},'Interpreter','none')
+            end
             % Annotation Box + save figure
             % ----------------------
-            if mod(j1,npan)==0 || (mod(j1,npan)~=0 && i==length(vars_irf))
+            if mod(j1,npan)==0 || (mod(j1,npan)~=0 && var_iter==length(vars_irf))
                 annotation('textbox', [0.1,0,0.35,0.05],'String', 'Linear','Color','Blue','horizontalalignment','center','interpreter','none');
                 annotation('textbox', [0.55,0,0.35,0.05],'String', 'Piecewise','Color','Red','horizontalalignment','center','interpreter','none');
-                dyn_saveas(hfig,[DirectoryName,filesep,M_.fname,'_irf_occbin_',simulname_,shocknames{j},'_',shocksign,'_',int2str(ifig)],options_.nodisplay,options_.graph_format);
+                dyn_saveas(hfig,[DirectoryName,filesep,M_.fname,'_irf_occbin_',save_name,shock_names{shock_iter},'_',shocksign,'_',int2str(ifig)],options_.nodisplay,options_.graph_format);
+                if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
+                    fprintf(fidTeX,'\\begin{figure}[H]\n');
+                    fprintf(fidTeX,'\\centering \n');
+                    fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_irf_occbin_%s}\n',options_.figures.textwidth*min((j1-1)/ncols,1),...
+                        [DirectoryName '/' ,M_.fname],[save_name,shock_names{shock_iter},'_',shocksign,'_',int2str(ifig)]);
+                    fprintf(fidTeX,'\\caption{OccBin IRFs to  %s shock to %s}\n',plot_title_sign,shock_names{shock_iter});
+                    fprintf(fidTeX,'\\label{Fig:occbin_irfs:%s}\n',[save_name,shock_names{shock_iter},'_',shocksign,'_',int2str(ifig)]);
+                    fprintf(fidTeX,'\\end{figure}\n');
+                    fprintf(fidTeX,'\n');
+                end
             end
-        end 
+        end
     end
 end
 
-warning('on','all')
+if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
+    fprintf(fidTeX,'%% End Of TeX file.');
+    fclose(fidTeX);
 end

matlab/+occbin/plot_regimes.m

@@ -1,4 +1,26 @@
 function plot_regimes(regimes,M_,options_)
+% plot_regimes(regimes,M_,options_)
+% Inputs: 
+% - regimes     	[structure]     OccBin regime information
+% - M_              [structure]     Matlab's structure describing the model
+% - options_        [structure]     Matlab's structure containing the options
+
+% Copyright © 2021-2023 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 nperiods = size(regimes,2);
 nconstr = length(fieldnames(regimes(1)))/2;
@@ -13,9 +35,15 @@ else
 end    
 
 GraphDirectoryName = CheckPath('graphs',M_.dname);
+fhandle = dyn_figure(options_.nodisplay,'Name',[M_.fname ': OccBin regimes']);
 
-fhandle = dyn_figure(options_.nodisplay,'Name',[M_.fname ' occbin regimes']);
+latexFolder = CheckPath('latex',M_.dname);
-
+if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
+    fidTeX = fopen([latexFolder '/' M_.fname '_occbin_regimes.tex'],'w');
+    fprintf(fidTeX,'%% TeX eps-loader file generated by occbin.plot_regimes.m (Dynare).\n');
+    fprintf(fidTeX,['%% ' datestr(now,0) '\n']);
+    fprintf(fidTeX,' \n');
+end
 
 for k=1:nconstr
     subplot(nconstr,1,k)
@@ -36,12 +64,23 @@ for k=1:nconstr
             end
         end
     end
-    title(['regime ' int2str(k)])
+    xlim([1 nperiods])
-    xlabel('historic period')
+    title(['Regime ' int2str(k)])
-    ylabel('regime expected start')
+    xlabel('Historic period')
+    ylabel('Regime: expected start')
 end
-annotation('textbox',[.25,0,.15,.05],'String','Unbinding','Color','blue');
+annotation('textbox',[.25,0,.15,.05],'String','Slack','Color','blue');
-annotation('textbox',[.65,0,.15,.05],'String','Binding','Color','red');
+annotation('textbox',[.65,0,.2,.05],'String','Binding','Color','red');
-
 
 dyn_saveas(fhandle,[GraphDirectoryName, filesep, M_.fname '_occbin_regimes'],options_.nodisplay,options_.graph_format);
+
+if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
+    fprintf(fidTeX,'\\begin{figure}[H]\n');
+    fprintf(fidTeX,'\\centering \n');
+    fprintf(fidTeX,'\\includegraphics[width=%2.2f\\textwidth]{%s_occbin_regimes}\n',options_.figures.textwidth,[GraphDirectoryName '/' M_.fname]);
+    fprintf(fidTeX,'\\caption{OccBin: regime evolution over time.}\n');
+    fprintf(fidTeX,'\\label{Fig:occbin_regimes}\n');
+    fprintf(fidTeX,'\\end{figure}\n');
+    fprintf(fidTeX,'\n');
+    fclose(fidTeX);
+end

matlab/+occbin/set_default_options.m

@@ -55,9 +55,7 @@ if ismember(flag,{'forecast','all'})
     options_occbin_.forecast.maxit=30;
     options_occbin_.forecast.qmc=0;
     options_occbin_.forecast.replic=0;
-    options_occbin_.forecast.sepath=0;    
     options_occbin_.forecast.SHOCKS0=[];
-    options_occbin_.forecast.treepath=1; % number of branches
 end
 
 if ismember(flag,{'irf','all'})    
@@ -98,17 +96,12 @@ end
 
 if ismember(flag,{'plot_irf','all'})
     options_occbin_.plot_irf.add_steadystate = 0;
-    options_occbin_.plot_irf.exo_names = [];
-    options_occbin_.plot_irf.endo_names = M_.endo_names;
-    options_occbin_.plot_irf.endo_names_long = [];
     options_occbin_.plot_irf.endo_scaling_factor = [];
     options_occbin_.plot_irf.grid            = true;
     options_occbin_.plot_irf.ncols            = 3;
     options_occbin_.plot_irf.nrows            = 3;
     options_occbin_.plot_irf.shocksigns = {'pos','neg'}; 
     options_occbin_.plot_irf.simulname='';
-    options_occbin_.plot_irf.threshold = 10^-6;
-    options_occbin_.plot_irf.tplot = [];
 end
 
 if ismember(flag,{'plot_shock_decomp','all'})

matlab/collapse_figures_in_tabgroup.m

@@ -0,0 +1,42 @@
+function collapse_figures_in_tabgroup
+
+% Copyright © 2023 Eduard Benet Cerda
+% Copyright © 2024 Dynare Team
+%
+% This file is part of Dynare.
+%
+% Dynare is free software: you can redistribute it and/or modify
+% it under the terms of the GNU General Public License as published by
+% the Free Software Foundation, either version 3 of the License, or
+% (at your option) any later version.
+%
+% Dynare is distributed in the hope that it will be useful,
+% but WITHOUT ANY WARRANTY; without even the implied warranty of
+% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+% GNU General Public License for more details.
+%
+% You should have received a copy of the GNU General Public License
+% along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
+
+% Create a new figure with results
+fig = uifigure(Name = 'Dynare Results');
+
+% Add a grid layout to make sure it spans the entire width
+g = uigridlayout(fig, [1,1], Padding = 0);
+
+% Add a tabgroup
+tg = uitabgroup(g);
+
+% Find all figures with Dynare Tag 
+f = findobj('-regexp','tag','dynare-figure');
+
+% Loop over all figures and reparent them to a tab. Avoid legends, they are
+% automatically tied.
+for j = 1 : numel(f)
+    t = uitab(tg);
+    types = arrayfun(@class, f(j).Children, 'UniformOutput', false);
+    idx = ismember(types, 'matlab.graphics.illustration.Legend'); % no need to reparent legends
+    set(f(j).Children(~idx),'Parent',t)
+    t.Title = f(j).Name;
+    delete(f(j))
+end

matlab/dyn_figure.m

@@ -12,7 +12,7 @@ function h = dyn_figure(nodisplay, varargin)
 % SPECIAL REQUIREMENTS
 %    none
 
-% Copyright © 2012-2017 Dynare Team
+% Copyright © 2012-2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -30,7 +30,7 @@ function h = dyn_figure(nodisplay, varargin)
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 if nodisplay
-    h = figure(varargin{:},'visible','off');
+    h = figure(varargin{:},'visible','off','tag','dynare-figure');
 else
-    h = figure(varargin{:});
+    h = figure(varargin{:},'tag','dynare-figure');
 end

matlab/dynare.m

@@ -122,7 +122,7 @@ if isempty(dot_location)
     fnamelength = length(fname);
     fname1 = [fname '.dyn'];
     d = dir(fname1);
-    if length(d) == 0
+    if isempty(d)
         fname1 = [fname '.mod'];
     end
     fname = fname1;
@@ -138,7 +138,7 @@ if fnamelength + length('.set_auxiliary_variables') > namelengthmax()
     error('Dynare: the name of your .mod file is too long, please shorten it')
 end
 
-if ~isempty(strfind(fname,filesep))
+if contains(fname,filesep)
     fprintf('\nIt seems you are trying to call a .mod file not located in the "Current Folder". This is not possible (the %s symbol is not allowed in the name of the .mod file).\n', filesep)
     [pathtomodfile,basename] = fileparts(fname);
     if exist(pathtomodfile,'dir')
@@ -297,17 +297,13 @@ if status
     error('Dynare: preprocessing failed')
 end
 
-if ~ isempty(find(abs(fname) == 46))
-    fname = fname(:,1:find(abs(fname) == 46)-1) ;
-end
-
 % We need to clear the driver (and only the driver, because the "clear all"
 % within the driver will clean the rest)
-clear(['+' fname '/driver'])
+clear(['+' fname(1:end-4) '/driver'])
 
 try
     cTic = tic;
-    evalin('base',[fname '.driver']);
+    evalin('base',[fname(1:end-4) '.driver']);
     cToc = toc(cTic);
     if nargout
         DynareInfo.time.compute = cToc;
@@ -315,7 +311,7 @@ try
 catch ME
     W = evalin('caller','whos');
     diary off
-    if ismember(fname,{W(:).name})
+    if ismember(fname(1:end-4),{W(:).name})
         error('Your workspace already contains a variable with the same name as the mod-file. You need to delete it or rename the mod-file.')
     else
         rethrow(ME)

matlab/dynare_graph.m

@@ -1,54 +0,0 @@
-function dyn_graph=dynare_graph(dyn_graph,y,tit,x)
-% function dynare_graph(y,tit,x)
-% graphs
-%
-% INPUT
-%   figure_name: name of the figures
-%   colors: line colors
-%
-% OUTPUT
-%   dyn_graph: structure with figure information
-%
-% SPECIAL REQUIREMENT
-%   none
-
-% Copyright © 2006-2017 Dynare Team
-%
-% This file is part of Dynare.
-%
-% Dynare is free software: you can redistribute it and/or modify
-% it under the terms of the GNU General Public License as published by
-% the Free Software Foundation, either version 3 of the License, or
-% (at your option) any later version.
-%
-% Dynare is distributed in the hope that it will be useful,
-% but WITHOUT ANY WARRANTY; without even the implied warranty of
-% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-% GNU General Public License for more details.
-%
-% You should have received a copy of the GNU General Public License
-% along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
-
-if nargin < 4
-    x = (1:size(y,1))';
-end
-nplot = dyn_graph.plot_nbr + 1;
-if nplot > dyn_graph.max_nplot
-    figure('Name',dyn_graph.figure_name);
-    nplot = 1;
-end
-dyn_graph.plot_nbr = nplot;
-subplot(dyn_graph.nr,dyn_graph.nc,nplot);
-
-line_types = dyn_graph.line_types;
-line_type = line_types{1};
-for i=1:size(y,2)
-    if length(line_types) > 1
-        line_type = line_types{i};
-    end
-
-    plot(x,y(:,i),line_type);
-    hold on
-end
-title(tit);
-hold off

matlab/dyntable.m

@@ -63,7 +63,7 @@ end
 value_format  = sprintf('%%%d.%df', val_width, val_precis);
 header_string_format  = sprintf('%%%ds', val_width);
 
-if length(title) > 0
+if ~isempty(title)
     fprintf('\n\n%s\n',title);
 end
 
@@ -72,7 +72,7 @@ if nargin==9
     disp(optional_header)
 end
 
-if length(headers) > 0
+if ~isempty(headers)
     hh_tbl = sprintf(label_format_leftbound , headers{1});
     for i=2:length(headers)
         hh_tbl  = [hh_tbl sprintf(header_string_format, headers{i})];

matlab/estimation/dsge_likelihood.m

@@ -1,6 +1,6 @@
 function [fval,info,exit_flag,DLIK,Hess,SteadyState,trend_coeff,M_,options_,bayestopt_,dr] = dsge_likelihood(xparam1,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,BoundsInfo,dr, endo_steady_state, exo_steady_state, exo_det_steady_state,derivatives_info)
 % [fval,info,exit_flag,DLIK,Hess,SteadyState,trend_coeff,M_,options_,bayestopt_,oo_] = dsge_likelihood(xparam1,dataset_,dataset_info,options_,M_,estim_params_,bayestopt_,BoundsInfo,oo_,derivatives_info)
-% Evaluates the posterior kernel of a DSGE model using the specified
+% Evaluates the negative of the posterior kernel of a DSGE model using the specified
 % kalman_algo; the resulting posterior includes the 2*pi constant of the
 % likelihood function
 %
@@ -21,7 +21,7 @@ function [fval,info,exit_flag,DLIK,Hess,SteadyState,trend_coeff,M_,options_,baye
 % - derivatives_info    [structure]     derivative info for identification
 %
 % OUTPUTS
-% - fval                    [double]        scalar, value of the likelihood or posterior kernel.
+% - fval                    [double]        scalar, value of minus the likelihood or posterior kernel.
 % - info                    [integer]       4×1 vector, informations resolution of the model and evaluation of the likelihood.
 % - exit_flag               [integer]       scalar, equal to 1 (no issues when evaluating the likelihood) or 0 (not able to evaluate the likelihood).
 % - DLIK                    [double]        Vector with score of the likelihood
@@ -37,7 +37,7 @@ function [fval,info,exit_flag,DLIK,Hess,SteadyState,trend_coeff,M_,options_,baye
 % This function calls: dynare_resolve, lyapunov_symm, lyapunov_solver, compute_Pinf_Pstar, kalman_filter_d, missing_observations_kalman_filter_d,
 % univariate_kalman_filter_d, kalman_steady_state, get_perturbation_params_deriv, kalman_filter, missing_observations_kalman_filter, univariate_kalman_filter, priordens
 
-% Copyright © 2004-2023 Dynare Team
+% Copyright © 2004-2024 Dynare Team
 %
 % This file is part of Dynare.
 %

matlab/forecast_graphs.m

@@ -76,7 +76,7 @@ for j= 1:nvar
             fprintf(fidTeX,' \n');
         end
         n_fig =n_fig+1;
-        eval(['hh_fig=dyn_figure(options_.nodisplay,''Name'',''Forecasts (' int2str(n_fig) ')'');']);
+        hh_fig=dyn_figure(options_.nodisplay,'Name',['Forecasts (' int2str(n_fig) ')']);
         m = 1;
     end
     subplot(nr, nc, m);
@@ -138,7 +138,7 @@ if isfield(oo_.forecast,'HPDinf_ME')
                 fprintf(fidTeX,' \n');
             end
             n_fig =n_fig+1;
-            eval(['hh_fig=dyn_figure(options_.nodisplay,''Name'',''Forecasts (' int2str(n_fig) ')'');']);
+            hh_fig=dyn_figure(options_.nodisplay,'Name',['Forecasts including ME (' int2str(n_fig) ')']);
             m = 1;
         end
         subplot(nr,nc,m);

matlab/latex/collect_latex_files.m

@@ -35,7 +35,7 @@ fprintf(fid,'%s \n','\usepackage{psfrag}');
 fprintf(fid,'%s \n','\usepackage{graphicx}');
 fprintf(fid,'%s \n','\usepackage{epstopdf}');
 fprintf(fid,'%s \n','\usepackage{longtable,booktabs}');
-fprintf(fid,'%s \n','\usepackage{amsmath,amsfonts}');
+fprintf(fid,'%s \n','\usepackage{amsmath,amsfonts,amssymb}');
 fprintf(fid,'%s \n','\usepackage{breqn}');
 fprintf(fid,'%s \n','\usepackage{float,morefloats,caption}');
 fprintf(fid,'%s \n','\begin{document}');

matlab/model_info.m

@@ -280,20 +280,10 @@ function print_line(names,var_index,lead_lag,M_)
     else
         aux_index=find([M_.aux_vars(:).endo_index]==var_index);
         aux_type=M_.aux_vars(aux_index).type;
-        if ~isfield(M_.aux_vars(aux_index),'orig_lead_lag') || isempty(M_.aux_vars(aux_index).orig_lead_lag)
+        if lead_lag==0
-            if ismember(aux_type,[1,3])
+            str = subst_auxvar(var_index, [], M_);
-                str = subst_auxvar(var_index, -1, M_);
-            elseif ismember(aux_type,[0,2])
-                str = subst_auxvar(var_index, 1, M_);
-            else
-                if lead_lag==0
-                    str = subst_auxvar(var_index, [], M_);
-                else
-                    str = subst_auxvar(var_index, lead_lag, M_);
-                end
-            end
         else
-            str = subst_auxvar(var_index, M_.aux_vars(aux_index).orig_lead_lag, M_);
+            str = subst_auxvar(var_index, lead_lag, M_);
         end
         aux_orig_expression=M_.aux_vars(aux_index).orig_expr;
         if isempty(aux_orig_expression)

matlab/moments/compute_variance_decomposition.m

@@ -56,7 +56,7 @@ else
         var_decomp(stationary_vars,i) = vx2;
         variance_sum_loop = variance_sum_loop +vx2; %track overall variance over shocks
     end
-    if ~options_.pruning && max(abs(variance_sum_loop-var_stationary)./var_stationary) > 1e-4
+    if ~options_.pruning && max(abs(variance_sum_loop-var_stationary)./var_stationary) > 1e-4 && max(abs(variance_sum_loop-var_stationary))>1e-7
         warning(['Aggregate variance and sum of variances by shocks ' ...
             'differ by more than 0.01 %'])
     end

matlab/pac-tools/hVectors.m

@@ -18,7 +18,7 @@ function [h, lrcp] = hVectors(params, H, auxmodel, kind, id)
 %    params(2:end-1) ⟶ Autoregressive parameters.
 %    params(end)     ⟶ Discount factor.
 
-% Copyright © 2018-2021 Dynare Team
+% Copyright © 2018-2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -52,21 +52,21 @@ n = length(H);
 tmp = eye(n*m)-kron(G, transpose(H)); % inv(W2)
 
 switch kind
-  case 'll'
+  case 'll' % (A.84), page 28 in Brayton, Davis and Tulip (2000)
     h = A_1*A_b*((kron(iota(m, m), H))'*(tmp\kron(iota(m, m), iota(n, id))));
-  case 'dd'
+  case 'dd' % (A.79), page 26 in Brayton, Davis and Tulip (2000)
     h = A_1*A_b*(kron(iota(m, m)'*inv(eye(m)-G), H')*(tmp\kron(iota(m, m), iota(n, id))));
-  case 'dl'
+  case 'dl' % (A.74), page 24 in Brayton, Davis and Tulip (2000)
     h = A_1*A_b*(kron(iota(m, m)'*inv(eye(m)-G), (H'-eye(length(H))))*(tmp\kron(iota(m, m), iota(n, id))));
   otherwise
     error('Unknown kind value in PAC model.')
 end
 
-if nargin>1
+if nargout>1
     if isequal(kind, 'll')
         lrcp = NaN;
     else
         d = A_1*A_b*(iota(m, m)'*inv((eye(m)-G)*(eye(m)-G))*iota(m, m));
         lrcp = (1-sum(params(2:end-1))-d);
     end
-end
+end

matlab/perfect-foresight-models/solve_block_decomposed_problem.m

@@ -14,7 +14,7 @@ function [y, success, maxerror, per_block_status] = solve_block_decomposed_probl
 %   maxerror         [double]    ∞-norm of the residual
 %   per_block_status [struct]    vector structure with per-block information about convergence
 
-% Copyright © 2020-2023 Dynare Team
+% Copyright © 2020-2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -33,18 +33,21 @@ function [y, success, maxerror, per_block_status] = solve_block_decomposed_probl
 
 cutoff = 1e-15;
 
-if options_.stack_solve_algo==0
+switch options_.stack_solve_algo
-    mthd='Sparse LU';
+    case 0
-elseif options_.stack_solve_algo==1 || options_.stack_solve_algo==6
+        mthd='Sparse LU on stacked system';
-    mthd='LBJ';
+    case {1,6}
-elseif options_.stack_solve_algo==2
+        mthd='LBJ with LU solver';
-    mthd='GMRES';
+    case 2
-elseif options_.stack_solve_algo==3
+        mthd='GMRES on stacked system';
-    mthd='BICGSTAB';
+    case 3
-elseif options_.stack_solve_algo==4
+        mthd='BiCGStab on stacked system';
-    mthd='OPTIMPATH';
+    case 4
-else
+        mthd='Sparse LU solver with optimal path length on stacked system';
-    mthd='UNKNOWN';
+    case 7
+        mthd='Solver from solve_algo option on stacked system'
+    otherwise
+        error('Unsupported stack_solve_algo value')
 end
 if options_.verbosity
     printline(41)

matlab/send_endogenous_variables_to_workspace.m

@@ -20,5 +20,5 @@ function send_endogenous_variables_to_workspace()
 global M_ oo_
 
 for idx = 1:M_.endo_nbr
-    assignin('base', M_.endo_names{idx}, oo_.endo_simul(idx,:))
+    assignin('base', M_.endo_names{idx}, oo_.endo_simul(idx,:)')
 end

matlab/send_exogenous_variables_to_workspace.m

@@ -21,5 +21,5 @@ function send_exogenous_variables_to_workspace()
 global M_ oo_
 
 for idx = 1:M_.exo_nbr
-    assignin('base', M_.exo_names{idx}, oo_.exo_simul(:,idx))
+    assignin('base', M_.exo_names{idx}, oo_.exo_simul(:,idx)')
 end

matlab/stochastic_solver/k_order_pert.m

@@ -1,7 +1,7 @@
 function [dr,info] = k_order_pert(dr,M_,options_)
 % Compute decision rules using the k-order DLL from Dynare++
 
-% Copyright © 2009-2023 Dynare Team
+% Copyright © 2009-2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -30,6 +30,16 @@ if M_.maximum_endo_lead == 0 && order>1
            'backward models'])
 end
 
+if options_.aim_solver
+    error('Option aim_solver is not compatible with k_order_solver')
+end
+if options_.dr_cycle_reduction
+    error('Option dr=cycle_reduction is not compatible with k_order_solver')
+end
+if options_.dr_logarithmic_reduction
+    error('Option dr=logarithmic_reduction is not compatible with k_order_solver')
+end
+
 try
     [dynpp_derivs, dyn_derivs] = k_order_perturbation(dr,M_,options_);
 catch ME

matlab/utilities/dataset/makedataset.m

@@ -242,8 +242,14 @@ if newdatainterface
     end
 else
     % ...  or check that nobs is smaller than the number of observations in dataset_.
-    if nobs>dataset_.nobs
+    if FIRSTOBS>dataset_.dates(1)
-        error('makedataset: nobs (%s) cannot be greater than the last date in the dataset (%s)!', num2str(nobs), num2str(dataset_.nobs))
+        if FIRSTOBS+nobs-1>dataset_.dates(end)
+            error('makedataset: Given first_obs=%u and %u total observations in the dataset, the current nobs of %s must not be greater than %s!', options_.first_obs, dataset_.nobs, num2str(nobs), num2str(dataset_.nobs-find(dataset_.dates==FIRSTOBS)+1))
+        end
+    else
+        if nobs>dataset_.nobs
+            error('makedataset: nobs (%s) cannot be greater than the last date in the dataset (%s)!', num2str(nobs), num2str(dataset_.nobs))
+        end
     end
 end
 

matlab/utilities/general/display_parameter_values.m

@@ -1,17 +1,6 @@
-function dynare_graph_close()
+function display_parameter_values
-% function dynare_graph_close()
-% close a figure
-%
-% INPUT
-%   none
-%
-% OUTPUT
-%   none
-%
-% SPECIAL REQUIREMENT
-%   none
 
-% Copyright © 2006-2017 Dynare Team
+% Copyright © 2024 Dynare Team
 %
 % This file is part of Dynare.
 %
@@ -27,3 +16,10 @@ function dynare_graph_close()
 %
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
+global M_
+my_title='Current parameter values:';
+labels = M_.param_names;
+headers = {'Parameter'; 'Value'};
+lh = cellofchararraymaxlength(labels)+2;
+options_.noprint=false;
+dyntable(options_, my_title, headers, labels, M_.params, lh, 10, 6);

meson.build

@@ -21,9 +21,14 @@ c_compiler = meson.get_compiler('c')
 
 subdir('preprocessor/src')
 
-# Compatibility symlink
+# NB: the following two symlinks my be improved if this wishlist item is done:
-# NB: the following my be improved if this wishlist item is done:
 #  https://github.com/mesonbuild/meson/issues/11519
+
+install_symlink('dynare-preprocessor', install_dir : 'lib/dynare/preprocessor/',
+                pointing_to : '../../../bin/dynare-preprocessor'
+                              + (host_machine.system() == 'windows' ? '.exe' : ''))
+
+# Compatibility symlink
 install_symlink('dynare_m' + (host_machine.system() == 'windows' ? '.exe' : ''),
                 install_dir : 'lib/dynare/matlab/preprocessor64',
                 pointing_to : '../../../../bin/dynare-preprocessor'
@@ -34,17 +39,12 @@ install_symlink('dynare_m' + (host_machine.system() == 'windows' ? '.exe' : ''),
 install_subdir('matlab', install_dir : 'lib/dynare',
                exclude_files : [ 'utilities/tests/.git' ,
                                  'utilities/tests/.gitignore',
-                                 'modules/reporting/.git',
+                                 'dseries/.git',
-                                 'modules/reporting/.gitignore',
+                                 'dseries/.gitignore',
-                                 'modules/reporting/.gitlab-ci.yml',
+                                 'dseries/src/modules/matlab-fame-io/.git',
-                                 'modules/dseries/.git',
+                                 'dseries/src/modules/matlab-fame-io/.gitignore',
-                                 'modules/dseries/.gitignore',
+                                 'dseries/.gitmodules',
-                                 'modules/dseries/src/modules/matlab-fame-io/.git',
+                                 'dseries/.gitlab-ci.yml' ])
-                                 'modules/dseries/src/modules/matlab-fame-io/.gitignore',
-                                 'modules/dseries/src/modules/matlab-fame-io/doc/.gitignore',
-                                 'modules/dseries/src/modules/matlab-fame-io/tests/FameDatabases/.gitignore',
-                                 'modules/dseries/.gitmodules',
-                                 'modules/dseries/.gitlab-ci.yml' ])
 
 sed_exe = find_program('sed')
 custom_target(output : 'dynare_version.m', input : 'matlab/dynare_version.m.in',
@@ -53,6 +53,9 @@ custom_target(output : 'dynare_version.m', input : 'matlab/dynare_version.m.in',
               install : true,
               install_dir : 'lib/dynare/matlab')
 
+install_subdir('contrib/ms-sbvar/TZcode/MatlabFiles',
+               install_dir : 'lib/dynare/contrib/ms-sbvar/TZcode')
+
 ### MEX files
 
 mex_incdir = include_directories('mex/sources')
@@ -157,7 +160,20 @@ if get_option('build_for') == 'matlab'
   umfpack_dep = declare_dependency(link_args : '-lmwumfpack', dependencies : blas_dep)
   ut_dep = declare_dependency(link_args : '-lut')
 
-  slicot_dep = declare_dependency(dependencies : [ fortran_compiler.find_library('slicot64_pic'), blas_dep, lapack_dep ])
+  # Workaround for Meson bug https://github.com/mesonbuild/meson/issues/12757
+  # Use the C compiler as a fallback for detecting SLICOT under Linux with
+  # prefer_static=true (but still try the Fortran compiler to honour the -B
+  # option in fortran_args, as documented). Needed for building the MATLAB
+  # Online package.
+  if get_option('prefer_static') and host_machine.system() == 'linux'
+    slicot_dep_tmp = fortran_compiler.find_library('slicot64_pic', required : false)
+    if not slicot_dep_tmp.found()
+      slicot_dep_tmp = c_compiler.find_library('slicot64_pic')
+    endif
+    slicot_dep = declare_dependency(dependencies : [ slicot_dep_tmp, blas_dep, lapack_dep ])
+  else
+    slicot_dep = declare_dependency(dependencies : [ fortran_compiler.find_library('slicot64_pic'), blas_dep, lapack_dep ])
+  endif
 else # Octave build
   octave_exe = find_program('octave', required : not meson.is_cross_build(), disabler : true)
   mkoctfile_exe = find_program('mkoctfile')
@@ -291,8 +307,11 @@ if get_option('prefer_static')
 
   # NB: constructing a dependency object with link_args : ['-Wl,-Bstatic', '-lgomp', '-Wl,-Bdynamic'] does not work,
   # because it reorders the three arguments and puts -lgomp at the end
-  openmp_dep_tmp = cpp_compiler.find_library('gomp', static : true)
+  if host_machine.system() != 'linux'
-  openmp_dep = declare_dependency(dependencies : [ openmp_dep, openmp_dep_tmp ])
+    # Under Debian 12, trying to link (static) libgomp.a in a MEX fails.
+    openmp_dep_tmp = cpp_compiler.find_library('gomp', static : true)
+    openmp_dep = declare_dependency(dependencies : [ openmp_dep, openmp_dep_tmp ])
+  endif
 endif
 
 # For use when creating intermediate static libraries to be incorporated in MEX files
@@ -1892,7 +1911,9 @@ endforeach
 git_exe = find_program('git', required : false)
 etags_exe = find_program('etags', required : false)
 
-if git_exe.found() and etags_exe.found()
+fs = import('fs')
+
+if fs.is_dir('.git') and git_exe.found() and etags_exe.found()
   all_files = run_command(git_exe,
                           [ '--git-dir=@0@/.git'.format(meson.project_source_root()),
                             'ls-files', '--recurse-submodules',

mex/sources/bytecode/Interpreter.cc

@@ -1,5 +1,5 @@
 /*
- * Copyright © 2007-2023 Dynare Team
+ * Copyright © 2007-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -4710,23 +4710,26 @@ Interpreter::Simulate_Newton_Two_Boundaries(
           switch (stack_solve_algo)
             {
             case 0:
-              mexPrintf("MODEL SIMULATION: (method=Sparse LU)\n");
+              mexPrintf("MODEL SIMULATION: (method=Sparse LU solver on stacked system)\n");
               break;
             case 2:
-              mexPrintf(preconditioner_print_out("MODEL SIMULATION: (method=GMRES)\n",
+              mexPrintf(
-                                                 preconditioner, false)
+                  preconditioner_print_out("MODEL SIMULATION: (method=GMRES on stacked system)\n",
-                            .c_str());
+                                           preconditioner, false)
+                      .c_str());
               break;
             case 3:
-              mexPrintf(preconditioner_print_out("MODEL SIMULATION: (method=BiCGStab)\n",
+              mexPrintf(preconditioner_print_out(
-                                                 preconditioner, false)
+                            "MODEL SIMULATION: (method=BiCGStab on stacked system)\n",
+                            preconditioner, false)
                             .c_str());
               break;
             case 4:
-              mexPrintf("MODEL SIMULATION: (method=Sparse LU & optimal path length)\n");
+              mexPrintf("MODEL SIMULATION: (method=Sparse LU solver with optimal path length on "
+                        "stacked system)\n");
               break;
             case 5:
-              mexPrintf("MODEL SIMULATION: (method=Sparse Gaussian Elimination)\n");
+              mexPrintf("MODEL SIMULATION: (method=LBJ with Sparse Gaussian Elimination)\n");
               break;
             }
         }

mex/sources/libkorder/sylv/QuasiTriangular.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004-2011 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -240,16 +240,11 @@ public:
     real = r;
   }
   virtual ~_matrix_iter() = default;
-  bool
+  [[nodiscard]] bool
   operator==(const _Self& it) const
   {
     return ptr == it.ptr;
   }
-  bool
-  operator!=(const _Self& it) const
-  {
-    return ptr != it.ptr;
-  }
   _TRef
   operator*() const
   {

mex/sources/libkorder/sylv/Vector.cc

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004-2011 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -114,42 +114,6 @@ Vector::Vector(mxArray* p) :
     throw SYLV_MES_EXCEPTION("This is not a dense array of real doubles.");
 }
 
-bool
-Vector::operator==(const Vector& y) const
-{
-  return ConstVector(*this) == y;
-}
-
-bool
-Vector::operator!=(const Vector& y) const
-{
-  return ConstVector(*this) != y;
-}
-
-bool
-Vector::operator<(const Vector& y) const
-{
-  return ConstVector(*this) < y;
-}
-
-bool
-Vector::operator<=(const Vector& y) const
-{
-  return ConstVector(*this) <= y;
-}
-
-bool
-Vector::operator>(const Vector& y) const
-{
-  return ConstVector(*this) > y;
-}
-
-bool
-Vector::operator>=(const Vector& y) const
-{
-  return ConstVector(*this) >= y;
-}
-
 void
 Vector::zeros()
 {
@@ -323,17 +287,10 @@ ConstVector::operator==(const ConstVector& y) const
   return i == len;
 }
 
-bool
+std::partial_ordering
-ConstVector::operator<(const ConstVector& y) const
+ConstVector::operator<=>(const ConstVector& y) const
 {
-  int i = std::min(len, y.len);
+  return std::lexicographical_compare_three_way(data, data + len, y.data, y.data + y.len);
-  int ii = 0;
-  while (ii < i && operator[](ii) == y[ii])
-    ii++;
-  if (ii < i)
-    return operator[](ii) < y[ii];
-  else
-    return len < y.len;
 }
 
 double

mex/sources/libkorder/sylv/Vector.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004-2011 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -25,6 +25,7 @@
    to avoid running virtual method invokation mechanism. Some
    members, and methods are thus duplicated */
 
+#include <compare>
 #include <complex>
 #include <utility>
 
@@ -108,15 +109,6 @@ public:
     return s;
   }
 
-  // Exact equality.
-  bool operator==(const Vector& y) const;
-  bool operator!=(const Vector& y) const;
-  // Lexicographic ordering.
-  bool operator<(const Vector& y) const;
-  bool operator<=(const Vector& y) const;
-  bool operator>(const Vector& y) const;
-  bool operator>=(const Vector& y) const;
-
   virtual ~Vector()
   {
     if (destroy)
@@ -219,29 +211,9 @@ public:
     return s;
   }
   // Exact equality
-  bool operator==(const ConstVector& y) const;
+  [[nodiscard]] bool operator==(const ConstVector& y) const;
-  bool
-  operator!=(const ConstVector& y) const
-  {
-    return !operator==(y);
-  }
   // Lexicographic ordering
-  bool operator<(const ConstVector& y) const;
+  [[nodiscard]] std::partial_ordering operator<=>(const ConstVector& y) const;
-  bool
-  operator<=(const ConstVector& y) const
-  {
-    return operator<(y) || operator==(y);
-  }
-  bool
-  operator>(const ConstVector& y) const
-  {
-    return !operator<=(y);
-  }
-  bool
-  operator>=(const ConstVector& y) const
-  {
-    return !operator<(y);
-  }
 
   [[nodiscard]] double getNorm() const;
   [[nodiscard]] double getMax() const;

mex/sources/libkorder/tl/equivalence.cc

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -35,12 +35,14 @@ OrdSequence::operator[](int i) const
    orderings can be used for different problem sizes. We order them
    according to the average, and then according to the first item. */
 
-bool
+std::partial_ordering
-OrdSequence::operator<(const OrdSequence& s) const
+OrdSequence::operator<=>(const OrdSequence& s) const
 {
   double ta = average();
   double sa = s.average();
-  return (ta < sa || ((ta == sa) && (operator[](0) > s[0])));
+  if (auto cmp1 = ta <=> sa; cmp1 != 0)
+    return cmp1;
+  return operator[](0) <=> s[0];
 }
 
 bool

mex/sources/libkorder/tl/equivalence.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -55,6 +55,7 @@
 
 #include "int_sequence.hh"
 
+#include <compare>
 #include <list>
 #include <string>
 #include <vector>
@@ -62,7 +63,7 @@
 /* Here is the abstraction for an equivalence class. We implement it as
    vector<int>. We have a constructor for empty class, copy
    constructor. What is important here is the ordering operator
-   operator<() and methods for addition of an integer, and addition of
+   operator<=>() and methods for addition of an integer, and addition of
    another sequence. Also we provide method has() which returns true if a
    given integer is contained. */
 
@@ -74,9 +75,9 @@ public:
   OrdSequence() : data()
   {
   }
-  bool operator==(const OrdSequence& s) const;
+  [[nodiscard]] bool operator==(const OrdSequence& s) const;
   int operator[](int i) const;
-  bool operator<(const OrdSequence& s) const;
+  [[nodiscard]] std::partial_ordering operator<=>(const OrdSequence& s) const;
   [[nodiscard]] const std::vector<int>&
   getData() const
   {
@@ -120,12 +121,7 @@ public:
   // Copy constructor plus gluing i1 and i2 in one class
   Equivalence(const Equivalence& e, int i1, int i2);
 
-  bool operator==(const Equivalence& e) const;
+  [[nodiscard]] bool operator==(const Equivalence& e) const;
-  bool
-  operator!=(const Equivalence& e) const
-  {
-    return !operator==(e);
-  }
   [[nodiscard]] int
   getN() const
   {

mex/sources/libkorder/tl/gs_tensor.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -76,16 +76,11 @@ public:
   // Constructs the tensor dimensions for slicing (see the implementation for details)
   TensorDimens(const IntSequence& ss, const IntSequence& coor);
 
-  bool
+  [[nodiscard]] bool
   operator==(const TensorDimens& td) const
   {
     return nvs == td.nvs && sym == td.sym;
   }
-  bool
-  operator!=(const TensorDimens& td) const
-  {
-    return !operator==(td);
-  }
 
   [[nodiscard]] int
   dimen() const

mex/sources/libkorder/tl/int_sequence.cc

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -103,10 +103,10 @@ IntSequence::operator==(const IntSequence& s) const
   return std::equal(data, data + length, s.data, s.data + s.length);
 }
 
-bool
+std::strong_ordering
-IntSequence::operator<(const IntSequence& s) const
+IntSequence::operator<=>(const IntSequence& s) const
 {
-  return std::lexicographical_compare(data, data + length, s.data, s.data + s.length);
+  return std::lexicographical_compare_three_way(data, data + length, s.data, s.data + s.length);
 }
 
 bool

mex/sources/libkorder/tl/int_sequence.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -43,6 +43,7 @@
 #define INT_SEQUENCE_HH
 
 #include <algorithm>
+#include <compare>
 #include <initializer_list>
 #include <utility>
 #include <vector>
@@ -119,12 +120,7 @@ public:
     if (destroy)
       delete[] data;
   }
-  bool operator==(const IntSequence& s) const;
+  [[nodiscard]] bool operator==(const IntSequence& s) const;
-  bool
-  operator!=(const IntSequence& s) const
-  {
-    return !operator==(s);
-  }
   int&
   operator[](int i)
   {
@@ -141,15 +137,10 @@ public:
     return length;
   }
 
-  /* We provide two orderings. The first operator<() is the linear
+  /* We provide two orderings. The first operator<=>() is the linear
      lexicographic ordering, the second less() is the non-linear Cartesian
      ordering. */
-  bool operator<(const IntSequence& s) const;
+  [[nodiscard]] std::strong_ordering operator<=>(const IntSequence& s) const;
-  bool
-  operator<=(const IntSequence& s) const
-  {
-    return (operator==(s) || operator<(s));
-  }
   [[nodiscard]] bool lessEq(const IntSequence& s) const;
   [[nodiscard]] bool less(const IntSequence& s) const;
 

mex/sources/libkorder/tl/kron_prod.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -86,7 +86,7 @@ public:
 
   KronProdDimens& operator=(const KronProdDimens& kd) = default;
   KronProdDimens& operator=(KronProdDimens&& kd) = default;
-  bool
+  [[nodiscard]] bool
   operator==(const KronProdDimens& kd) const
   {
     return rows == kd.rows && cols == kd.cols;

mex/sources/libkorder/tl/permutation.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -95,7 +95,7 @@ public:
   Permutation(const Permutation& p, int i) : permap(p.permap.insert(p.size(), i))
   {
   }
-  bool
+  [[nodiscard]] bool
   operator==(const Permutation& p) const
   {
     return permap == p.permap;

mex/sources/libkorder/tl/ps_tensor.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -115,7 +115,7 @@ public:
   {
     per.apply(nvmax);
   }
-  bool
+  [[nodiscard]] bool
   operator==(const PerTensorDimens& td) const
   {
     return TensorDimens::operator==(td) && per == td.per;

mex/sources/libkorder/tl/symmetry.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -126,16 +126,11 @@ public:
   {
     return run;
   }
-  bool
+  [[nodiscard]] bool
-  operator=(const symiterator& it)
+  operator==(const symiterator& it)
   {
     return dim == it.dim && run == it.run;
   }
-  bool
-  operator!=(const symiterator& it)
-  {
-    return !operator=(it);
-  }
 };
 
 /* The class SymmetrySet defines a set of symmetries of the given length

mex/sources/libkorder/tl/tensor.hh

@@ -1,6 +1,6 @@
 /*
  * Copyright © 2004 Ondra Kamenik
- * Copyright © 2019-2023 Dynare Team
+ * Copyright © 2019-2024 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -150,16 +150,11 @@ public:
     {
       return offset;
     }
-    bool
+    [[nodiscard]] bool
     operator==(const index& n) const
     {
       return offset == n.offset;
     }
-    bool
-    operator!=(const index& n) const
-    {
-      return offset != n.offset;
-    }
     [[nodiscard]] const IntSequence&
     getCoor() const
     {

scripts/matlab-online/packageDynare.sh

@@ -4,12 +4,8 @@ set -exo pipefail
 # Creates a dynare-X.Y.mltbx in the current repository, using the settings below.
 # Needs to be run from Ubuntu 22.04 LTS, with the needed packages installed.
 
-DYNAREVER=5.5
 X13ASVER=1-1-b60
-LIBLOCATION="/MATLAB Add-Ons/Toolboxes/Dynare/mex/matlab/libs"
 MATLABPATH=/opt/MATLAB/R2023b
-# MatIO has been recompiled by hand to avoid the dependency on HDF5, which is a nightmare
-MATIO_PREFIX=/home/sebastien/usr
 # TODO: change size and put white background for better rendering in MATLAB Add-Ons browser
 DYNARE_PNG_LOGO=../../preprocessor/doc/logos/dlogo.png
 
@@ -21,48 +17,31 @@ cleanup ()
 }
 trap cleanup EXIT
 
-pushd "$tmpdir"
+pushd ../..
+meson setup -Dmatlab_path="$MATLABPATH" -Dbuildtype=release -Dprefer_static=true "$tmpdir"/build-matlab-online
 
-# Get Dynare
+cd "$tmpdir"/build-matlab-online
-wget -q https://www.dynare.org/release/source/dynare-$DYNAREVER.tar.xz
+meson compile
-tar xf dynare-$DYNAREVER.tar.xz
+meson install --destdir "$tmpdir"
+DYNAREVER=$(meson introspect --projectinfo | jq -r '.version')
 
-# Build Dynare
+cd ..
-cd dynare-$DYNAREVER
+strip usr/local/bin/dynare-preprocessor
-
+strip usr/local/lib/dynare/mex/matlab/*.mexa64
-# Use static libstdc++ otherwise the one shipped with MATLAB creates problem (since it overrides the systemd-wide one from Ubuntu)
-./configure --with-matlab="$MATLABPATH" --with-matio="$MATIO_PREFIX" --disable-octave --disable-doc --disable-dynare++ LDFLAGS=-static-libstdc++
-make -j$(nproc)
-
-strip preprocessor/dynare-preprocessor
-strip mex/matlab/*.mexa64
-
-# Patch mex files to look into ./lib folder
-for f in mex/matlab/*.mexa64
-do
-    patchelf --set-rpath "$LIBLOCATION" $f
-done
-
-# Grab the shared libraries needed
-mkdir mex/matlab/libs
-for l in libgsl.so.27 libgslcblas.so.0
-do
-    cp /usr/lib/x86_64-linux-gnu/$l mex/matlab/libs/
-    # Patch rpath to find dependencies at runtime
-    patchelf --set-rpath "$LIBLOCATION" mex/matlab/libs/$l
-done
-cp "$MATIO_PREFIX"/lib/libmatio.so.11 mex/matlab/libs/
-patchelf --set-rpath "$LIBLOCATION" mex/matlab/libs/libmatio.so.11 # Probably not needed
 
 # Get X13 binary from the Census Bureau website
 # The binary from Ubuntu has some shared library dependencies, so it is safer to use a static binary
 wget -q https://www2.census.gov/software/x-13arima-seats/x13as/unix-linux/program-archives/x13as_ascii-v${X13ASVER}.tar.gz
 tar xf x13as_ascii-v${X13ASVER}.tar.gz
-mkdir -p matlab/modules/dseries/externals/x13/linux/64
+
-cp x13as/x13as_ascii matlab/modules/dseries/externals/x13/linux/64/x13as
+# Populate staging area for the zip
+cp -pRL usr/local/lib/dynare dynare # -L is needed to dereference the preprocessor symlink
+mkdir -p dynare/matlab/dseries/externals/x13/linux/64
+cp -p x13as/x13as_ascii dynare/matlab/dseries/externals/x13/linux/64/x13as
 
 # zip dynare
-zip -q -r "$tmpdir"/dynare.zip contrib/jsonlab contrib/ms-sbvar/TZcode/MatlabFiles examples matlab mex/matlab preprocessor/dynare-preprocessor license.txt NEWS.md README.md VERSION
+cd dynare
+zip -q -r "$tmpdir"/dynare.zip *
 
 # make toolbox
 popd

tests/deterministic_simulations/purely_backward/ar1.mod

@@ -33,6 +33,6 @@ if ~oo_.deterministic_simulation.status
 end
 send_endogenous_variables_to_workspace;
 
-if max(abs(y'-[1; exp(cumprod([1; rho*ones(9, 1)]))]))>options_.dynatol.x
+if max(abs(y-[1; exp(cumprod([1; rho*ones(9, 1)]))]))>options_.dynatol.x
     error('Wrong solution!')
 end

tests/moments/example1_hp_test.mod

@@ -77,12 +77,12 @@ send_endogenous_variables_to_workspace;
 options_.nomoments=0;
 oo_unfiltered_all_shocks=oo_;
 
-[junk, y_filtered]=sample_hp_filter(y',1600);
+[junk, y_filtered]=sample_hp_filter(y,1600);
-[junk, c_filtered]=sample_hp_filter(c',1600);
+[junk, c_filtered]=sample_hp_filter(c,1600);
-[junk, k_filtered]=sample_hp_filter(k',1600);
+[junk, k_filtered]=sample_hp_filter(k,1600);
-[junk, a_filtered]=sample_hp_filter(a',1600);
+[junk, a_filtered]=sample_hp_filter(a,1600);
-[junk, h_filtered]=sample_hp_filter(h',1600);
+[junk, h_filtered]=sample_hp_filter(h,1600);
-[junk, b_filtered]=sample_hp_filter(b',1600);
+[junk, b_filtered]=sample_hp_filter(b,1600);
 
 verbatim;
 total_std_all_shocks_filtered_sim=std([y_filtered c_filtered k_filtered a_filtered h_filtered b_filtered]);
@@ -112,12 +112,12 @@ stoch_simul(order=1,nofunctions,hp_filter=0,periods=2500000,nomoments);
 send_endogenous_variables_to_workspace;
 oo_unfiltered_one_shock=oo_;
 
-[junk, y_filtered]=sample_hp_filter(y',1600);
+[junk, y_filtered]=sample_hp_filter(y,1600);
-[junk, c_filtered]=sample_hp_filter(c',1600);
+[junk, c_filtered]=sample_hp_filter(c,1600);
-[junk, k_filtered]=sample_hp_filter(k',1600);
+[junk, k_filtered]=sample_hp_filter(k,1600);
-[junk, a_filtered]=sample_hp_filter(a',1600);
+[junk, a_filtered]=sample_hp_filter(a,1600);
-[junk, h_filtered]=sample_hp_filter(h',1600);
+[junk, h_filtered]=sample_hp_filter(h,1600);
-[junk, b_filtered]=sample_hp_filter(b',1600);
+[junk, b_filtered]=sample_hp_filter(b,1600);
 
 verbatim;
 total_std_one_shock_filtered_sim=std([y_filtered c_filtered k_filtered a_filtered h_filtered b_filtered]);

tests/occbin/filter/NKM.mod

@@ -13,66 +13,66 @@
 
 // ----------------- Defintions -----------------------------------------//
 var        
-    c          //1  Consumption 
+    c       $c$         (long_name='Consumption')
-    n          //2  Labor
+    n                   (long_name='Labor')
-    y          //5  Output
+    y       $y$         (long_name='Output')
-    yf         //6  Final goods       
+    yf                  (long_name='Final goods')
-    yg         //11 Output growth gap
+    yg                  (long_name='Output growth gap')
-    w          //12 Real wage rate
+    w                   (long_name='Real wage rate')
-    wf         //13 Flexible real wage
+    wf                  (long_name='Flexible real wage')
-    pigap      //15 Inflation rate -> pi(t)/pibar = pigap
+    pigap               (long_name='Inflation rate -> pi(t)/pibar = pigap')
-    inom    ${i^{nom}}$       //16 Nominal interest rate
+    inom    ${i^{nom}}$ (long_name='Nominal interest rate')
-    inomnot    //17 Notional interest rate
+    inomnot             (long_name='Notional interest rate')
-    mc         //19 Real marginal cost
+    mc                  (long_name='Real marginal cost')
-    lam  ${\lambda}$      //20 Inverse marginal utility of wealth  
+    lam  ${\lambda}$    (long_name='Inverse marginal utility of wealth')
-    g          //21 Growth shock       
+    g                   (long_name='Growth shock')
-    s          //22 Risk premium shock
+    s                   (long_name='Risk premium shock')
-    mp         //23 Monetary policy shock    
+    mp                  (long_name='Monetary policy shock')
-    pi   ${\pi}$   //24 Observed inflation    
+    pi   ${\pi}$        (long_name='Observed inflation')
     @#if !(small_model)
-        x          //3  Investment
+        x               (long_name='Investment')
-        k          //4  Capital    
+        k               (long_name='Capital')
-        u          //7  Utilization cost
+        u               (long_name='Utilization cost')
-        ups        //8  Utilization choice
+        ups             (long_name='Utilization choice')
-        wg         //9  Real wage growth gap    
+        wg              (long_name='Real wage growth gap')    
-        xg         //10 Investment growth
+        xg              (long_name='Investment growth')
-        rk         //14 Real rental rate
+        rk              (long_name='Real rental rate')
-        q          //18 Tobins q   
+        q               (long_name='Tobins q')
     @#endif
 ;
 varexo          
-    epsg    ${\varepsilon_g}$   // Productivity growth shock
+    epsg    ${\varepsilon_g}$   (long_name='Productivity growth shock')
-    epsi       // Notional interest rate shock
+    epsi                        (long_name='Notional interest rate shock')
-    epss       // Risk premium shock
+    epss                        (long_name='Risk premium shock')
 ;        
 parameters
     // Calibrated Parameters    
-    beta    $\beta$    // Discount factor
+    beta    $\beta$    (long_name='Discount factor')
-    chi         // Labor disutility scale
+    chi                 (long_name='Labor disutility scale')
-    thetap      // Elasticity of subs. between intermediate goods
+    thetap              (long_name='Elasticity of subs. between intermediate goods')
-    thetaw      // Elasticity of subs. between labor types
+    thetaw              (long_name='Elasticity of subs. between labor types')
-    nbar        // Steady state labor
+    nbar                (long_name='Steady state labor')
-    eta         // Inverse frish elasticity of labor supply
+    eta                 (long_name='Inverse frish elasticity of labor supply')
-    delta       // Depreciation
+    delta               (long_name='Depreciation')
-    alpha       // Capital share
+    alpha               (long_name='Capital share')
-    gbar        // Mean growth rate
+    gbar                (long_name='Mean growth rate')
-    pibar       // Inflation target
+    pibar               (long_name='Inflation target')
-    inombar     // Steady gross nom interest rate    
+    inombar             (long_name='Steady gross nom interest rate')
-    inomlb      // Effective lower bound on gross nominal interest rate    
+    inomlb              (long_name='Effective lower bound on gross nominal interest rate')
-    sbar        // Average risk premium
+    sbar                (long_name='Average risk premium')
     // Parameters for DGP and Estimated parameters
-    varphip     // Rotemberg price adjustment cost
+    varphip             (long_name='Rotemberg price adjustment cost')
-    varphiw     // Rotemberg wage adjustment cost
+    varphiw             (long_name='Rotemberg wage adjustment cost')
-    h          	// Habit persistence
+    h          	        (long_name='Habit persistence')
-    rhos        // Persistence
+    rhos                (long_name='Persistence')
-    rhoi    	// Persistence
+    rhoi    	        (long_name='Persistence')
-    sigz        // Standard deviation technology
+    sigz                (long_name='Standard deviation technology')
-    sigs        // Standard deviation risk premia
+    sigs                (long_name='Standard deviation risk premia')
-    sigi        // Standard deviation mon pol
+    sigi                (long_name='Standard deviation mon pol')
-    phipi       // Inflation responsiveness
+    phipi               (long_name='Inflation responsiveness')
-    phiy        // Output responsiveness
+    phiy                (long_name='Output responsiveness')
-    nu          // Investment adjustment cost
+    nu                  (long_name='Investment adjustment cost')
-    sigups      // Utilization    
+    sigups              (long_name='Utilization')
  ;
 
 
@@ -316,36 +316,34 @@ varobs yg inom pi;
 
     // forecast starting from period 42, zero shocks (default)
     smoother2histval(period=42);
-    [oo, error_flag] = occbin.forecast(options_,M_,oo_,8);
+    [forecast, error_flag] = occbin.forecast(options_,M_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state,8);
     // forecast with stochastic shocks
     options_.occbin.forecast.qmc=true;
     options_.occbin.forecast.replic=127;
-    [oo1, error_flag] = occbin.forecast(options_,M_,oo_,8);
+    [forecast1, error_flag] = occbin.forecast(options_,M_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state,8);
+
+    figure('Name','OccBin: Forecasts')
+    subplot(2,1,1)
+    plot(1:8,forecast.piecewise.Mean.inom,'b-',1:8,forecast1.piecewise.Mean.inom,'r--')
+    subplot(2,1,2)
+    plot(1:8,forecast.piecewise.Mean.y,'b-',1:8,forecast1.piecewise.Mean.y,'r--')
 
     // GIRF given states in 42 and shocks in 43
     t0=42;
     options_.occbin.irf.exo_names=M_.exo_names;
     options_.occbin.irf.t0=t0;
-    oo_ = occbin.irf(M_,oo_,options_);
+    oo_.occbin.irfs = occbin.irf(M_,oo_,options_);
 
-    vars_irf = {
+    var_list_ = {'c','n','y','pigap','inom','inomnot'};
-    'c', 'consumption'    
-    'n', 'labor'    
-    'y', 'output'    
-    'pigap', 'inflation rate'
-    'inom', 'interest rate'  
-    'inomnot', 'shadow rate'
-    };
 
-    options_.occbin.plot_irf.exo_names = M_.exo_names;
-    options_.occbin.plot_irf.endo_names = vars_irf(:,1);
-    options_.occbin.plot_irf.endo_names_long = vars_irf(:,2);
     // if you want to scale ...
     // options_occbin_.plot_irf.endo_scaling_factor = vars_irf(:,3);
     options_.occbin.plot_irf.simulname = ['t0_' int2str(t0)];
-    options_.occbin.plot_irf.tplot = min(40,options_.irf);
+    options_.irf=40;
-    occbin.plot_irfs(M_,oo_,options_);
+    occbin.plot_irfs(M_,oo_.occbin.irfs,options_,var_list_);
-
+    var_list_={};
+    options_.occbin.plot_irf.simulname = ['t0_' int2str(t0) '_full'];
+    occbin.plot_irfs(M_,oo_.occbin.irfs,options_,var_list_);
     oo0=oo_;
     // use smoother_redux
     estimation(
@@ -374,7 +372,7 @@ varobs yg inom pi;
             consider_all_endogenous,heteroskedastic_filter,filter_step_ahead=[1],smoothed_state_uncertainty);
             
     // show initial condition effect of IF
-    figure,
+    figure('Name','OccBin: Smoothed shocks')
     subplot(221)
     plot([oo0.SmoothedShocks.epsg oo_.SmoothedShocks.epsg]), title('epsg')
     subplot(222)
@@ -382,7 +380,7 @@ varobs yg inom pi;
     subplot(223)
     plot([oo0.SmoothedShocks.epss oo_.SmoothedShocks.epss]), title('epss')
     legend('PKF','IF')
-    figure,
+    figure('Name','OccBin: Smoothed Variables')
     subplot(221)
     plot([oo0.SmoothedVariables.inom oo_.SmoothedVariables.inom]), title('inom')
     subplot(222)

windows/deps/Makefile

@@ -157,6 +157,7 @@ mingw64: tarballs/mingw-w64-x86_64-gcc-$(MINGW64_GCC_VERSION)-any.pkg.tar.zst ta
 	touch $@
 
 tarballs/mingw-w64-x86_64-%-any.pkg.tar.zst:
+	mkdir -p tarballs
 	wget $(WGET_OPTIONS) -O $@ http://repo.msys2.org/mingw/x86_64/$(notdir $@)
 
 clean-msys2:

windows/dynare.nsi

@@ -52,7 +52,7 @@ Section "Dynare core (preprocessor and M-files)"
  File README.txt ..\NEWS.md ..\license.txt
 
  SetOutPath $INSTDIR\matlab
- File /r ..\matlab\*.m ..\build-win-matlab\dynare_version.m
+ File /r ..\matlab\*.m ..\matlab\*.json ..\build-win-matlab\dynare_version.m
 
  SetOutPath $INSTDIR\preprocessor
  File ..\build-win-matlab\preprocessor\src\dynare-preprocessor.exe