
x_start=[0.01, 0.01, 0.01, 0.01]';
%define parameters to be optimized and their upper and lower bound
x_opt_name={'alpha_le',0,Inf
            'alpha_qky',0,Inf
            'alpha_ve',0,Inf
            'alpha_cg',0,Inf
            };

options_.nofunctions=1;
options_.nograph=1;
options_.verbosity=0;

%set noprint option to suppress error messages within optimizer
options_.noprint=1;
options_.TeX=0;
% set csminwel options
H0 = 1e-2*eye(length(x_start)); %Initial Hessian 
crit = 1e-8; %Tolerance
nit = 1000;  %Number of iterations
[fhat,x_opt_hat] = csminwel(@welfare_objective,x_start,H0,[],crit,nit,x_opt_name);



function [outvalue]=welfare_objective(x_opt,x_opt_name)
% function [outvalue]=welfare_objective(x_opt,x_opt_name)

global oo_ options_ M_

%% set parameter for use in Dynare
for ii=1:size(x_opt_name,1)
    set_param_value(x_opt_name{ii,1},x_opt(ii));
end

if any(x_opt<cell2mat(x_opt_name(:,2))) || any(x_opt>cell2mat(x_opt_name(:,3))) %make sure parameters are inside their bounds
    outvalue=10e6+sum([x_opt].^2); %penalty function
    return
end

var_list_ = char('omega_e');
info = stoch_simul(var_list_); %get decision rules and moments
if info(1) %filter out error code
    outvalue=1e5+sum([x_opt].^2);
    return;
end
outvalue=-oo_.mean(strmatch('omega_e',var_list_,'exact')); %extr