function plot_2D_afs(AFS, addfpip, D)
% plots an 2D AFS of a three component system with one, three or more 
% subsets
% 
% input: 
% AFS ........ area of feasible solutions
% 
% optional input:
% addfpip .... use "addfpip = 1" to plot also the outer and inner polygon,
% use "addfpip = 0" to skip them
% additional plot of the inner and the outer polygon requires data D
% 
% examples:
% plot_2D_afs(AFS)
% plot_2D_afs(newsegA,1,D)

col = [0     0.4470     0.7410;
    0.8500    0.3250    0.0980;
    0.9290    0.6940    0.1250;
    0.4940    0.1840    0.5560;
    0.4660    0.6740    0.1880;
    0.3010    0.7450    0.9330];

COL = [0.85     0.9     1;
    0.995    0.8250    0.75;
    1    0.9    0.75;
    0.85    0.6    0.925;
    0.75    0.95    0.6;
    0.7    0.9    1];

s = 3;

figure;
% computation of the inner polygon and the outer polygon if selected
if exist('addfpip') == 1 && exist('D') && addfpip == 1
    [U,S,V] = svds(D,s);
    for i=1:s
        if min(V(:,i)) < -max(V(:,i))
            V(:,i) = -V(:,i);
            U(:,i) = -U(:,i);
        end
    end
    a = diag(1./U(:,1)/S(1,1))*U(:,2:3)*S(2:3,2:3);
    b = diag(1./V(:,1))*V(:,2:3);
    cips = convhull(a);
    cipc = convhull(b);
    IPS = a(cips,:);
    IPC = b(cipc,:);
    FPS = [];
    m = size(IPC,1);
    for i=1:m-1
        X = [IPC(i,:); IPC(i+1,:)];
        FPS(end+1,:) = X\[-1;-1];
    end
    FPS(end+1,:) = FPS(1,:);
    
    plot(FPS(:,1), FPS(:,2), '--k', 'linewidth', 0.5)
    line(IPS(:,1), IPS(:,2), 'linestyle', '-.', 'color', 'k', 'linewidth', 0.5)
    hold on;
end

% plot the single subset(s)
if isempty(AFS{3})
    fill(AFS{1}(:,1), AFS{1}(:,2), [.7 .7 .7], 'linewidth', 1, 'edgecolor', 'k')
    hold on;
    fill(AFS{2}(:,1), AFS{2}(:,2), [1 1 1], 'linewidth', 1, 'edgecolor', 'k')
    if exist('IPS')
        line(IPS(:,1), IPS(:,2), 'linestyle', '-.', 'color', 'k', 'linewidth', 0.5)
    end
else
    for i=1:length(AFS)
        fill(AFS{i}(:,1), AFS{i}(:,2), COL(i,:), 'linewidth', 1, 'edgecolor', col(i,:));
        hold on;
    end
end
hold off
axis tight
xl = get(gca, 'xlim');
yl = get(gca, 'ylim');
xlabel('$x_1$','interpreter','latex','fontsize',12)
ylabel('$x_2$','interpreter','latex','fontsize',12)
title('Area of feasible solutions','FontWeight','Normal')
set(gca, 'xlim', 1.05*xl, 'ylim', 1.05*yl)
