
clc;  clear;  close
load ("data/Parameters.mat",'m');

JJ = 41;
dis = 0.025;

%% The equilibrium
% Case 1: competition + prohibit
tic
m.DataPolicy = 'prohibit';  % prohibit free
m.MarketStruct = 'competition';  % competition oligopoly monopoly 
m.Privacy = 'nonePrivacy';  % nonePrivacy privacyExternality

m.t_cu = 0;       % 国家U单边最优关税
m.tau_cu = m.d_cu*(1+m.t_cu);
for i = 1:JJ
    m.B_c = 0.5 + dis*(i-1);
    m.B_u = 1;
    ub = [];   lb = 0;   opt0 = 0;
    [opt_Case_uc, fval, exitflag] =  knitromatlab(@(opt)func_obj_t_uc(opt,m), opt0, [], [], [], [], lb, ub, [], [], [], 'ktropt_singlestart.opt');      
    m.t_uc = opt_Case_uc;
    m.tau_uc = m.d_uc*(1+m.t_uc);

    [wel1,Case11] = func_eqm_iter(m);   
    Case1_uc.t_uc(i,1) = m.t_uc;
    Case1_uc.wel_u(i,1) = Case11.wel_u;
end

m.t_uc = 0;       % 国家C单边最优关税
m.tau_uc = m.d_uc*(1+m.t_uc);
for i = 1:JJ
    m.B_c = 0.5 + dis*(i-1);
    m.B_u = 1;
    ub = [];   lb = 0;   opt0 = 0;
    [opt_Case_cu, fval, exitflag] =  knitromatlab(@(opt)func_obj_t_cu(opt,m), opt0, [], [], [], [], lb, ub, [], [], [], 'ktropt_singlestart.opt');      
    m.t_cu = opt_Case_cu;
    m.tau_cu = m.d_cu*(1+m.t_cu);

    [wel1,Case11] = func_eqm_iter(m);   
    Case1_cu.t_cu(i,1) = m.t_cu;
    Case1_cu.wel_c(i,1) = Case11.wel_c;
end


% Case 2: competition + free
m.DataPolicy = 'free';  % prohibit free
m.MarketStruct = 'competition';  % competition oligopoly monopoly 
m.Privacy = 'nonePrivacy';  % nonePrivacy privacyExternality

m.t_cu = 0;       % 国家U单边最优关税
m.tau_cu = m.d_cu*(1+m.t_cu);
for i = 1:JJ
    m.B_c = 0.5 + dis*(i-1);
    m.B_u = 1;
    ub = [];   lb = 0;   opt0 = 0;
    [opt_Case_uc, fval, exitflag] =  knitromatlab(@(opt)func_obj_t_uc(opt,m), opt0, [], [], [], [], lb, ub, [], [], [], 'ktropt_singlestart.opt');      
    m.t_uc = opt_Case_uc;
    m.tau_uc = m.d_uc*(1+m.t_uc);

    [wel2,Case22] = func_eqm_iter(m);   
    Case2_uc.t_uc(i,1) = m.t_uc;
    Case2_uc.wel_u(i,1) = Case22.wel_u;
end

m.t_uc = 0;       % 国家C单边最优关税
m.tau_uc = m.d_uc*(1+m.t_uc);
for i = 1:JJ
    m.B_c = 0.5 + dis*(i-1);
    m.B_u = 1;
    ub = [];   lb = 0;   opt0 = 0;
    [opt_Case_cu, fval, exitflag] =  knitromatlab(@(opt)func_obj_t_cu(opt,m), opt0, [], [], [], [], lb, ub, [], [], [], 'ktropt_singlestart.opt');      
    m.t_cu = opt_Case_cu;
    m.tau_cu = m.d_cu*(1+m.t_cu);

    [wel2,Case22] = func_eqm_iter(m);   
    Case2_cu.t_cu(i,1) = m.t_cu;
    Case2_cu.wel_c(i,1) = Case22.wel_c;
end

B_c = 0.5:dis:0.5+dis*(JJ-1);
B_c = B_c';

%% Figure：算力和福利
warning('off','all')

set(0, 'defaulttextinterpreter','latex')
set(0, 'defaultAxesTickLabelInterpreter','latex'); 
set(0, 'defaultLegendInterpreter','latex');
set(groot,'defaultLineLineWidth',2); 

TightPlot.ColumeNumber = 2; 
TightPlot.RowNumber = 2;
TightPlot.GapW = 0.1;
TightPlot.GapH = 0.15;
TightPlot.MarginsLower = 0.15;
TightPlot.MarginsUpper = 0.1;
TightPlot.MarginsLeft = 0.1;
TightPlot.MarginsRight = 0.1;

picture1 = figure('color',[1 1 1]);
p = tight_subplot(TightPlot.ColumeNumber,TightPlot.RowNumber,...
    [TightPlot.GapH TightPlot.GapW],...
    [TightPlot.MarginsLower TightPlot.MarginsUpper],...
    [TightPlot.MarginsLeft TightPlot.MarginsRight]); 

axes(p(1));
line(B_c,Case1_cu.t_cu, 'Color', 'k', 'LineStyle', '-', 'LineWidth', 1.25);
hold on;
scatter(B_c, Case1_uc.t_uc, 6, '^', 'k');
hold off;
set(gca, 'XTickLabelMode', 'auto');
set(gca,'yTickLabel',num2str(get(gca,'yTick')','%.2f'));
xlabel('$ B_{C} $','Interpreter','latex', 'FontSize',10);
ylabel('最优关税','interpreter','latex', 'FontSize',10);
title(p(1),'（a）竞争/禁止流动','FontSize',10)

axes(p(2));
line(B_c,Case2_cu.t_cu, 'Color', 'k', 'LineStyle', '-', 'LineWidth', 1.25);
hold on;
scatter(B_c, Case2_uc.t_uc, 6, '^', 'k');
hold off;
set(gca, 'XTickLabelMode', 'auto');
set(gca,'yTickLabel',num2str(get(gca,'yTick')','%.2f'));
xlabel('$ B_{C} $','Interpreter','latex', 'FontSize',10);
ylabel('最优关税','interpreter','latex', 'FontSize',10);
title(p(2),'（b）竞争/自由流动','FontSize',10)

axes(p(3));
line(B_c,Case1_cu.wel_c, 'Color', 'k', 'LineStyle', '-', 'LineWidth', 1.25);
hold on;
scatter(B_c, Case1_uc.wel_u, 6, '^', 'k');
hold off;
set(gca, 'XTickLabelMode', 'auto');
set(gca,'yTickLabel',num2str(get(gca,'yTick')','%.2f'));
xlabel('$ B_{C} $','Interpreter','latex', 'FontSize',10);
ylabel('福利','interpreter','latex', 'FontSize',10);
title(p(3),'（c）竞争/禁止流动','FontSize',10)

axes(p(4));
line(B_c,Case2_cu.wel_c, 'Color', 'k', 'LineStyle', '-', 'LineWidth', 1.25);
hold on;
scatter(B_c, Case2_uc.wel_u, 6, '^', 'k');
hold off;
set(gca, 'XTickLabelMode', 'auto');
set(gca,'yTickLabel',num2str(get(gca,'yTick')','%.2f'));
xlabel('$ B_{C} $','Interpreter','latex', 'FontSize',10);
ylabel('福利','interpreter','latex', 'FontSize',10);
title(p(4),'（d）竞争/自由流动','FontSize',10)

handle = legend({'国家 C', '国家 U','最优关税'});
set(handle,'Position', [0.45 0.0000001 0.1 0.05],'Orientation','horizon','Interpreter','latex','FontSize',8)

fprintf('\n=================================================================\n');
fprintf('         图III4 国家𝐶不同算力水平下的最优关税和相应福利（𝜌 = 0.5）\n');
fprintf('=================================================================\n\n');
print(picture1,'result/图III4.png','-r600','-djpeg');

