diff --git a/matlab/check_bounds_and_definiteness_estimation.m b/matlab/check_bounds_and_definiteness_estimation.m
new file mode 100644
index 0000000000000000000000000000000000000000..88b3f8460581204bb89f509fbb32502035410a32
--- /dev/null
+++ b/matlab/check_bounds_and_definiteness_estimation.m
@@ -0,0 +1,109 @@
+function [fval,info,exit_flag,M_,Q,H]=check_bounds_and_definiteness_estimation(xparam1, M_, options_, estim_params_, bayestopt_)
+% function [fval,info,exit_flag]=check_bounds_and_definiteness_estimation(xparam1, M_, options_, estim_params_, bayestopt_)
+% Checks whether parameter vector satisfies 
+%
+% INPUTS
+% - xparam1                 [double]              n by 1 vector, estimated parameters.
+% - M_                      [struct]              Matlab's structure describing the Model.
+% - options_                [struct]              Matlab's structure describing the options.
+% - estim_params_           [struct]              Matlab's structure describing the estimated_parameters.
+% - bayestopt_              [struct]              Matlab's structure specifying the bounds on the paramater values (initialized by dynare,aka bayesopt_).
+%
+% OUTPUTS
+% - fval                    [double]              scalar, value of the likelihood or posterior kernel.
+% - info                    [integer]             4 by 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).
+% - Q                       [matrix]              Covariance matrix of structural shocks
+% - H                       [matrix]              Covariance matrix of measurement errors
+
+% Copyright (C) 2020 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 <http://www.gnu.org/licenses/>.
+
+fval        = [];
+exit_flag   = 1;
+info        = zeros(4,1);
+Q=[];
+H=[];
+% Return, with endogenous penalty, if some parameters are smaller than the lower bound of the prior domain.
+if any(xparam1<bayestopt_.lb)
+    k = find(xparam1(:) < bayestopt_.lb);
+    fval = Inf;
+    exit_flag = 0;
+    info(1) = 41;
+    info(4) = sum((bayestopt_.lb(k)-xparam1(k)).^2);
+    return
+end
+
+% Return, with endogenous penalty, if some parameters are greater than the upper bound of the prior domain.
+if any(xparam1>bayestopt_.ub)
+    k = find(xparam1(:)>bayestopt_.ub);
+    fval = Inf;
+    exit_flag = 0;
+    info(1) = 42;
+    info(4) = sum((xparam1(k)-bayestopt_.ub(k)).^2);
+    return
+end
+
+M_ = set_all_parameters(xparam1,estim_params_,M_);
+
+Q = M_.Sigma_e;
+H = M_.H;
+
+if ~issquare(Q) || estim_params_.ncx || isfield(estim_params_,'calibrated_covariances')
+    [Q_is_positive_definite, penalty] = ispd(Q(estim_params_.Sigma_e_entries_to_check_for_positive_definiteness,estim_params_.Sigma_e_entries_to_check_for_positive_definiteness));
+    if ~Q_is_positive_definite
+        fval = Inf;
+        exit_flag = 0;
+        info(1) = 43;
+        info(4) = penalty;
+        return
+    end
+    if isfield(estim_params_,'calibrated_covariances')
+        correct_flag=check_consistency_covariances(Q);
+        if ~correct_flag
+            penalty = sum(Q(estim_params_.calibrated_covariances.position).^2);
+            fval = Inf;
+            exit_flag = 0;
+            info(1) = 71;
+            info(4) = penalty;
+            return
+        end
+    end
+
+end
+
+if ~issquare(H) || estim_params_.ncn || isfield(estim_params_,'calibrated_covariances_ME')
+    [H_is_positive_definite, penalty] = ispd(H(estim_params_.H_entries_to_check_for_positive_definiteness,estim_params_.H_entries_to_check_for_positive_definiteness));
+    if ~H_is_positive_definite
+        fval = Inf;
+        exit_flag = 0;
+        info(1) = 44;
+        info(4) = penalty;
+        return
+    end
+    if isfield(estim_params_,'calibrated_covariances_ME')
+        correct_flag=check_consistency_covariances(H);
+        if ~correct_flag
+            penalty = sum(H(estim_params_.calibrated_covariances_ME.position).^2);
+            fval = Inf;
+            exit_flag = 0;
+            info(1) = 72;
+            info(4) = penalty;
+            return
+        end
+    end
+end
\ No newline at end of file
diff --git a/matlab/dsge_likelihood.m b/matlab/dsge_likelihood.m
index 7b568277c201d5bc8f86b3fc41d531e896f803c4..d43988af6665a444f4f35b2330e9244d473132a1 100644
--- a/matlab/dsge_likelihood.m
+++ b/matlab/dsge_likelihood.m
@@ -181,85 +181,11 @@ end
 % 1. Get the structural parameters & define penalties
 %------------------------------------------------------------------------------
 
-% Return, with endogenous penalty, if some parameters are smaller than the lower bound of the prior domain.
-if isestimation(DynareOptions) && ~isequal(DynareOptions.mode_compute,1) && any(xparam1<BoundsInfo.lb)
-    k = find(xparam1<BoundsInfo.lb);
-    fval = Inf;
-    exit_flag = 0;
-    info(1) = 41;
-    info(4)= sum((BoundsInfo.lb(k)-xparam1(k)).^2);
-    if analytic_derivation
-        DLIK=ones(length(xparam1),1);
-    end
-    return
-end
-
-% Return, with endogenous penalty, if some parameters are greater than the upper bound of the prior domain.
-if isestimation(DynareOptions) && ~isequal(DynareOptions.mode_compute,1) && any(xparam1>BoundsInfo.ub)
-    k = find(xparam1>BoundsInfo.ub);
-    fval = Inf;
-    exit_flag = 0;
-    info(1) = 42;
-    info(4)= sum((xparam1(k)-BoundsInfo.ub(k)).^2);
-    if analytic_derivation
-        DLIK=ones(length(xparam1),1);
-    end
+[fval,info,exit_flag,Model,Q,H]=check_bounds_and_definiteness_estimation(xparam1, Model, DynareOptions, EstimatedParameters, BoundsInfo);
+if info(1)
     return
 end
 
-% Get the diagonal elements of the covariance matrices for the structural innovations (Q) and the measurement error (H).
-Model = set_all_parameters(xparam1,EstimatedParameters,Model);
-
-Q = Model.Sigma_e;
-H = Model.H;
-
-% Test if Q is positive definite.
-if ~issquare(Q) || EstimatedParameters.ncx || isfield(EstimatedParameters,'calibrated_covariances')
-    [Q_is_positive_definite, penalty] = ispd(Q(EstimatedParameters.Sigma_e_entries_to_check_for_positive_definiteness,EstimatedParameters.Sigma_e_entries_to_check_for_positive_definiteness));
-    if ~Q_is_positive_definite
-        fval = Inf;
-        exit_flag = 0;
-        info(1) = 43;
-        info(4) = penalty;
-        return
-    end
-    if isfield(EstimatedParameters,'calibrated_covariances')
-        correct_flag=check_consistency_covariances(Q);
-        if ~correct_flag
-            penalty = sum(Q(EstimatedParameters.calibrated_covariances.position).^2);
-            fval = Inf;
-            exit_flag = 0;
-            info(1) = 71;
-            info(4) = penalty;
-            return
-        end
-    end
-end
-
-% Test if H is positive definite.
-if ~issquare(H) || EstimatedParameters.ncn || isfield(EstimatedParameters,'calibrated_covariances_ME')
-    [H_is_positive_definite, penalty] = ispd(H(EstimatedParameters.H_entries_to_check_for_positive_definiteness,EstimatedParameters.H_entries_to_check_for_positive_definiteness));
-    if ~H_is_positive_definite
-        fval = Inf;
-        info(1) = 44;
-        info(4) = penalty;
-        exit_flag = 0;
-        return
-    end
-    if isfield(EstimatedParameters,'calibrated_covariances_ME')
-        correct_flag=check_consistency_covariances(H);
-        if ~correct_flag
-            penalty = sum(H(EstimatedParameters.calibrated_covariances_ME.position).^2);
-            fval = Inf;
-            exit_flag = 0;
-            info(1) = 72;
-            info(4) = penalty;
-            return
-        end
-    end
-end
-
-
 %------------------------------------------------------------------------------
 % 2. call model setup & reduction program
 %------------------------------------------------------------------------------
diff --git a/matlab/dsge_var_likelihood.m b/matlab/dsge_var_likelihood.m
index d68ed2149fa0fe7bf05cfb9a2816022d3f458825..69e64899c17515ea9938c477b69697b1590edf5f 100644
--- a/matlab/dsge_var_likelihood.m
+++ b/matlab/dsge_var_likelihood.m
@@ -101,38 +101,11 @@ mYX = evalin('base', 'mYX');
 mXY = evalin('base', 'mXY');
 mXX = evalin('base', 'mXX');
 
-% Return, with endogenous penalty, if some dsge-parameters are smaller than the lower bound of the prior domain.
-if isestimation(DynareOptions) && DynareOptions.mode_compute ~= 1 && any(xparam1 < BoundsInfo.lb)
-    k = find(xparam1 < BoundsInfo.lb);
-    fval = Inf;
-    exit_flag = 0;
-    info(1) = 41;
-    info(4)= sum((BoundsInfo.lb(k)-xparam1(k)).^2);
-    return
-end
-
-% Return, with endogenous penalty, if some dsge-parameters are greater than the upper bound of the prior domain.
-if isestimation(DynareOptions) && DynareOptions.mode_compute ~= 1 && any(xparam1 > BoundsInfo.ub)
-    k = find(xparam1 > BoundsInfo.ub);
-    fval = Inf;
-    exit_flag = 0;
-    info(1) = 42;
-    info(4) = sum((xparam1(k)-BoundsInfo.ub(k)).^2);
+[fval,info,exit_flag,Model,Q]=check_bounds_and_definiteness_estimation(xparam1, Model, DynareOptions, EstimatedParameters, BoundsInfo);
+if info(1)
     return
 end
 
-% Get the variance of each structural innovation.
-Q = Model.Sigma_e;
-for i=1:EstimatedParameters.nvx
-    k = EstimatedParameters.var_exo(i,1);
-    Q(k,k) = xparam1(i)*xparam1(i);
-end
-offset = EstimatedParameters.nvx;
-
-% Update Model.params and Model.Sigma_e.
-Model.params(EstimatedParameters.param_vals(:,1)) = xparam1(offset+1:end);
-Model.Sigma_e = Q;
-
 % Get the weight of the dsge prior.
 dsge_prior_weight = Model.params(dsge_prior_weight_idx);
 
diff --git a/matlab/non_linear_dsge_likelihood.m b/matlab/non_linear_dsge_likelihood.m
index 42d563a9bc58adc3a0c4717c7bad0adbc812cb89..1c4694af9ea9101f0ee15beed089e56e5a79fa19 100644
--- a/matlab/non_linear_dsge_likelihood.m
+++ b/matlab/non_linear_dsge_likelihood.m
@@ -63,77 +63,11 @@ end
 % 1. Get the structural parameters & define penalties
 %------------------------------------------------------------------------------
 
-% Return, with endogenous penalty, if some parameters are smaller than the lower bound of the prior domain.
-if isestimation(DynareOptions) && (DynareOptions.mode_compute~=1) && any(xparam1<BoundsInfo.lb)
-    k = find(xparam1(:) < BoundsInfo.lb);
-    fval = Inf;
-    exit_flag = 0;
-    info(1) = 41;
-    info(4) = sum((BoundsInfo.lb(k)-xparam1(k)).^2);
-    return
-end
-
-% Return, with endogenous penalty, if some parameters are greater than the upper bound of the prior domain.
-if isestimation(DynareOptions) && (DynareOptions.mode_compute~=1) && any(xparam1>BoundsInfo.ub)
-    k = find(xparam1(:)>BoundsInfo.ub);
-    fval = Inf;
-    exit_flag = 0;
-    info(1) = 42;
-    info(4) = sum((xparam1(k)-BoundsInfo.ub(k)).^2);
+[fval,info,exit_flag,Model,Q,H]=check_bounds_and_definiteness_estimation(xparam1, Model, DynareOptions, EstimatedParameters, BoundsInfo);
+if info(1)
     return
 end
 
-Model = set_all_parameters(xparam1,EstimatedParameters,Model);
-
-Q = Model.Sigma_e;
-H = Model.H;
-
-if ~issquare(Q) || EstimatedParameters.ncx || isfield(EstimatedParameters,'calibrated_covariances')
-    [Q_is_positive_definite, penalty] = ispd(Q(EstimatedParameters.Sigma_e_entries_to_check_for_positive_definiteness,EstimatedParameters.Sigma_e_entries_to_check_for_positive_definiteness));
-    if ~Q_is_positive_definite
-        fval = Inf;
-        exit_flag = 0;
-        info(1) = 43;
-        info(4) = penalty;
-        return
-    end
-    if isfield(EstimatedParameters,'calibrated_covariances')
-        correct_flag=check_consistency_covariances(Q);
-        if ~correct_flag
-            penalty = sum(Q(EstimatedParameters.calibrated_covariances.position).^2);
-            fval = Inf;
-            exit_flag = 0;
-            info(1) = 71;
-            info(4) = penalty;
-            return
-        end
-    end
-
-end
-
-if ~issquare(H) || EstimatedParameters.ncn || isfield(EstimatedParameters,'calibrated_covariances_ME')
-    [H_is_positive_definite, penalty] = ispd(H(EstimatedParameters.H_entries_to_check_for_positive_definiteness,EstimatedParameters.H_entries_to_check_for_positive_definiteness));
-    if ~H_is_positive_definite
-        fval = Inf;
-        exit_flag = 0;
-        info(1) = 44;
-        info(4) = penalty;
-        return
-    end
-    if isfield(EstimatedParameters,'calibrated_covariances_ME')
-        correct_flag=check_consistency_covariances(H);
-        if ~correct_flag
-            penalty = sum(H(EstimatedParameters.calibrated_covariances_ME.position).^2);
-            fval = Inf;
-            exit_flag = 0;
-            info(1) = 72;
-            info(4) = penalty;
-            return
-        end
-    end
-
-end
-
 %------------------------------------------------------------------------------
 % 2. call model setup & reduction program
 %------------------------------------------------------------------------------
diff --git a/matlab/utilities/general/isestimation.m b/matlab/utilities/general/isestimation.m
deleted file mode 100644
index 9e17d4b03a70a847067ba735e87f69f51df1f434..0000000000000000000000000000000000000000
--- a/matlab/utilities/general/isestimation.m
+++ /dev/null
@@ -1,26 +0,0 @@
-function b = isestimation(option)
-
-% Returns true if we are currently estimating a model.
-
-% Copyright (C) 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 <http://www.gnu.org/licenses/>.
-
-b = false;
-
-if ischar(option.mode_compute) || option.mode_compute || option.mh_replic
-    b = true;
-end
\ No newline at end of file