setpathdynare6_3
clear


if exist('param_overwrite.mat','file')
    delete param_overwrite.mat
end

if exist('param_overwrite_extra.mat','file')
    delete param_overwrite_extra.mat
end


mod_name = 'BBCDLoptrule';

eval(['dynare ',mod_name,' noclearall nowarn'])

% nash results 0.91  99.99
% cooperation results 0.95 99.99

param_list1 = char('c1gamrn','c1gamdpc');
param_init_val1 = [      0.95       99.99 ];

param_upper_bound1 = [0.999 100 ];
param_lower_bound1 = [0     1.1 ];


param_list2 = char('c2gamrn','c2gamdpc');
param_init_val2 = [     0.95          99.99]; 

param_upper_bound2 = [0.999  100];
param_lower_bound2 = [0      1.1];

objective_label = 'Welf';
objective_label1 = 'Welf1';
objective_label2 = 'Welf2';

shocks_labels = char('c1zerr','c2zerr','c2taupvarerr','c1taupvarerr');

tol0 = 0.001;       % tolerance for convergence of Nash loop
max_iter = 10;      % maximum number of iterations for the Nash loop


%% 

nparams1 = size(param_list1,1);
nparams2 = size(param_list2,1);

% test objective function
% objective_value = objective(param_init_val1,param_list1,param_upper_bound1,...
%         param_lower_bound1,mod_name,objective_label,shocks_labels);

    
% merge parameters and other choices for cooperative solution    
param_init_val = [param_init_val1, param_init_val2];
param_list = char(param_list1,param_list2);
param_upper_bound = [param_upper_bound1, param_upper_bound2];
param_lower_bound = [param_lower_bound1, param_lower_bound2];

options = optimset('Display','Iter','MaxFunEvals',1e4,'MaxIter',1e4,'TolFun',1e-5,'TolX',1e-3);
param_opt_cooperation=fminsearch('objective',param_init_val,options,param_list,param_upper_bound,param_lower_bound,mod_name,objective_label,shocks_labels);



%% set up the recursion for the Nash game
current_diff = 100;
this_iter = 0;

options = optimset('Display','Iter','MaxFunEvals',1e4,'MaxIter',200,'TolFun',1e-5,'TolX',1e-3);

while current_diff>tol0 && this_iter<=max_iter
    this_iter = this_iter+1;
    
    % optimize the parameters for the first rule keeping the parameters for
    % the second rule at their previously optimized values (or at their
    % initial values).
    
    % overwrite parameter values for rule 2
    for indx = 1:nparams2
        eval([param_list2(indx,:),'=param_init_val2(indx);'])
    end
    
    % save parameters for other rule to disk
    paramline = strmerge(cellstr(param_list2),' ');
    eval(['save param_overwrite_extra ',paramline])
    
    % optimize parameters for first rule
    param_opt_nash1=fminsearch('objective',param_init_val1,options,param_list1,param_upper_bound1,param_lower_bound1,mod_name,objective_label1,shocks_labels);
    
    % overwrite parameter values for rule 2
    for indx = 1:nparams1
        eval([param_list1(indx,:),'=param_init_val1(indx);'])
    end
    
    % save parameters for other rule to disk
    paramline = strmerge(cellstr(param_list1),' ');
    eval(['save param_overwrite_extra ',paramline])  
    
    % optimize the parameters for the other rule
    param_opt_nash2=fminsearch('objective',param_init_val2,options,param_list2,param_upper_bound2,param_lower_bound2,mod_name,objective_label2,shocks_labels);
    
    current_diff = max(max(abs(param_opt_nash1-param_init_val1)),...
                       max(abs(param_opt_nash2-param_init_val2)) );
    
    param_init_val1 = param_opt_nash1;
    param_init_val2 = param_opt_nash2;
end

%% clean up
delete param_overwrite
delete param_overwrite_extra