Plot array

clear;clc;

file='test';
q=0.9;
a=0.9;
c0=0.1;
c1=0.001;
X=2^(1/2);

for i=1:1:101
    for j=1:1:101

       b=0.01*i-0.01;
       D=(j-1)/100;

        p1 = (q*(a + b))/2 - (q*(q*(a + b) - 2*c0))/(8*(q - c1));
        p2 = q*(D + a) - (q*X*(q*(a + b) - 2*c0))/(8*(q - c1));
        w1 = (q*(2*c0 + a*q + b*q - 2*c1*(a + b)))/(4*(q - c1));
        w2f = -(q*(q*(a + b) - 2*c0 - 4*(q - c1)*(a + b)))/(8*(q - c1)) - (q*(q*(a + b) - 2*c0)^2)/(4*(q - c1)*(q*(a + b) - 2*c0 + 4*(q - c1)*(a - b + 2*D)));
        w2a = q*(D + a) - (X*q*(q*(a + b) - 2*c0))/(4*(q - c1));
        pmf = (q*(a + b))/2 - (q*(q*(a + b) - 2*c0 + 2*q*(q - b + 2*D)))/(8*(q - c1)) - (q^2*(q*(a + b)))/(8*(q - c1)^2);
        pma = q*(b - D) + (q*(c0 - q*(b - D)))/(2*(q - c1)) - (X*q*c1*(q*(a + b) - 2*c0))/(8*(q - c1)^2);
        c = c0 - c1*(a + b - 2*p1/q);
        c2f = c0 - c1*(a + b - p1/q - pmf/q);
        c2a = c0 - c1*(a + b - p2/q - pma/q);
        Po1 = (q*(a + b) - 2*c0)^2/(16*(q - c1));

        if a - b > 0 & q - c1 > 0 & q*(a + b) - 2*c0 > 0 & w1 - c > 0 & p1 - w1 > 0 & a - p1/q > 0 & b - p1/q > 0
            S1 = (q*(a + b) - 2*c0)^2/(16*(q - c1));
        else
            S1 = 0;
        end

        if a - b > 0 & q - c1 > 0 & q*(a + b) - 2*c0 > 0 & w1 - c > 0 & p1 - w1 > 0 & a - p1/q > 0 & b - p1/q > 0 & w2f - c2f > 0 & pmf - c2f > 0 & a + D - p1/q > 0 & b - D - pmf/q > 0 & p1 - w2f > 0
            F2 = Po1 + (q^2*(2*c0 + 4*D*c1 - 4*D*q + 2*a*c1 - 2*b*c1 - 3*a*q + b*q)^2)/(64*(q - c1)^3);
        else
            F2 = 0;
        end

        if a - b > 0 & q - c1 > 0 & q*(a + b) - 2*c0 > 0 & w1 - c > 0 & p1 - w1 > 0 & a - p1/q > 0 & b - p1/q > 0 & w2a - c2a > 0 & pma - c2a > 0 & a + D - p2/q > 0 & b - D - pma/q > 0 & p2 - w2a > 0
            A2 = - ((c0 + D*q - b*q)*(4*c0*c1 - 4*c0*q - 4*D*q^2 + 4*b*q^2 + 4*D*c1*q - 2*X*c0*c1 - 4*b*c1*q + X*a*c1*q + X*b*c1*q))/(16*(c1 - q)^2) - (X*((q^2*(4*b*c1 - 4*D*c1 - 8*c0 + X*a*c1 + X*b*c1) - 4*c0*c1^2 + 2*c1*q*(6*c0 + 2*D*c1 - X*c0 - 2*b*c1))/(8*(c1 - q)^2) + q*(D + a) + (X*q*(a*q - 2*c0 + b*q))/(4*(c1 - q)))*(a*q - 2*c0 + b*q))/(8*(c1 - q));
        else
            A2 = 0;
        end

        Y = [S1 F2 A2];
        Z = max(Y);
        if Z == 0

            %             fprintf('%4.2f,%4.2f,N/A\n',b,D)
            AA=[{'N/A'},NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,NaN,];
              str = {'N/A'};
        else
            if Z == S1
                %                 fprintf('%4.2f,%4.2f,S1\n',b,D)
                AA=[{'S1'},p1,NaN,w1,NaN,NaN,NaN,NaN,c,NaN,NaN,S1,NaN,NaN];
              str = {'S1'};
            end
            if Z == F2
                %                 fprintf('%4.2f,%4.2f,F2\n',b,D)
                AA=[{'F2'},p1,NaN,NaN,w2f,NaN,pmf,NaN,NaN,c2f,NaN,NaN,F2,NaN];
                  str = {'F2'};
            end
            if Z == A2
                %                 fprintf('%4.2f,%4.2f,A2\n',b,D)
                AA=[{'A2'},NaN,p2,NaN,NaN,w2a,NaN,pma,NaN,NaN,c2a,NaN,NaN,A2];
                str = {'A2'};
            end
        end
      matrix((i-1)*101+j,:)=[q,a,c0,c1,b,D,str,AA(2),AA(3),AA(4),AA(5),AA(6)...
          ,AA(7),AA(8),AA(9),AA(10),AA(11),AA(12),AA(13),AA(14)];

    end
end
title={'q','Alpha','c0','c1','Beta','Delta','Max','p1','p2','w1'...
    ,'w2f','w2a','pmf','pma','c','c2f','c2a','S1','F2','A2'};
sheet=1;
xlRange = 'A2';
xlswrite(file,title,sheet,'A1');
xlswrite(file,matrix,sheet,xlRange);