
clear all;
warning off;

%%%
% Custom parameters :
%
% degv : polynomials degrees to search
% family type : 1 - Fekete Points; 2 - Lebesgue Points
%%%

degv=2:5; 
family_type=2;

%%%
% Differential Evolution parameters :
% de_parms(1) is max number of iterations
% de_parms(2) is neighbourhood size
% de_parms(3) is value to reach
% de_parms(4) is gamma_F
% de_parms(5) is gamma_CR
%%%
de_parms=[150 1/50 1 0.5 0.9]; 

%%%
% Domain : 1 for SQUARE, 2 for DISK, 3 for INTERVAL
%%%
domain = 2;

%--------------------------------------------------------------------------
% PROGRAM FLOW STARTS HERE.
%--------------------------------------------------------------------------
s='SQUARE_LEB';
diaryname=cat(2,'diary_cb_',date,s,'.m');

diary on; diary(diaryname);

fprintf('function [pts, pts_orig,leb]=set_de_square_leb(deg)');
fprintf('\n \n \t switch deg \n ');

summary = [];

for iii=1:length(degv)


    deg=degv(iii);
	fprintf('\n case %d\n',deg);

    if length(family_type) > 0
        fprintf('\n \t %% STAGE 2');

		%%%
		% Domain-specific Initialization
		%%%
        [pts,V_leb] = init(deg,domain);

		%%%
		% Points transformation
		%%%
		x_in = f(pts,domain);

		%%%
		% Parameter setting
		%%%
        F_VTR=de_parms(3);
        I_D=length(x_in);
        F_weight = de_parms(4);
        F_CR = de_parms(5);
        mydelta=de_parms(2);
        xlow=(x_in-mydelta)';
        xup=(x_in+mydelta)';
        maxit=de_parms(1);

		%%%
		% Core DE routine
		%%%
        [x_new,S_y,I_nf]=myRundeopt(F_VTR,I_D,F_weight,F_CR,xlow,...
            xup,maxit,deg,V_leb,...
            family_type,domain);
        if size(x_new,1) < size(x_new,2)
            x_new = x_new';
        end
		pts_in = pts;
		pts_out = g(x_new, domain);

		%%%
		% Selection
		%%%
        value_in = target(pts_in,domain,V_leb,deg);
        value_out=target(pts_out,domain,V_leb,deg);

        fprintf('\n \t %% IN: %2.5e OUT: %2.5e DIFF: %1.1e',...
            value_in,value_out,abs(value_in-value_out));

        if value_in < value_out
            pts_new=pts_in; value_new=value_in;
        else
            pts_new=pts_out; value_new=value_out;
        end

        %%%%%%
        % Output Section
        %%%%%%

        fprintf('\n \t pts_orig=[');
        for ggg=1:length(pts)
			fprintf('\n');
			for hhh=1:size(pts,2)
                fprintf('\t%1.15e',pts(ggg,hhh));
			end
        end
        fprintf('\n \t ];');
        fprintf('\n \t pts=[');
        for ggg=1:length(pts_new)
			fprintf('\n');
			for hhh=1:size(pts,2)
                fprintf('\t%1.15e',pts_new(ggg,hhh));
			end
        end
        fprintf('\n \t ];');
		clf;
		if size(pts,2) <= 1
		    plot(pts,pts_new);
		else
		    hold on;
		    plot(pts(:,1),pts(:,2),"ro");
		    plot(pts_new(:,1),pts_new(:,2),"gx");
    	    hold off;
		end
        end

        fprintf('\n\t leb= %2.5e ;\n',value_new);

		%%%
		% Stats Section
		%%%
        row = size(summary)(1) +1;
		V_old = vand(deg,pts,domain);
		V_new = vand(deg,pts_new,domain);
		leb_old = norm(V_old'\V_leb',1);
		leb_new = norm(V_new'\V_leb',1);
		lagr = max(max(abs(V_new'\V_leb')));
		summary(row,:) = [deg,leb_old,leb_new,lagr,cond(V_new,1),abs(det(V_new))];

   pause;
end

fprintf('\n%%%%%%%%%%');
fprintf('\n%% Summary table');
fprintf('\n%% Degree:\tLeb. old\tLeb. new\tMax Lagr.\tCond.\tV Det.');
for i=1:size(summary)(1)
    fprintf('\n%% %d\t',summary(i,1));
	for j=2:size(summary)(2)
		fprintf('%2.5e\t',summary(i,j));
	end
end
fprintf('\n%%%%%%%%%%');


diary off;
