function table = run_MOPB09_ex4()
% function table = run_MOPB09_ex4()
% synthesize gain-scheduling state-feedback gains for the Example 4 presented in
% (Montagner, Oliveira, Peres & Bliman, IJC, 2009) using the conditions of Theorem 3 
% given in the same paper. The system is given originally in the affine form and was
% borrowed from (de Souza & Trofino, IJRNC, 2006). For each relaxation level (degree g
% and Polya's level d) a bissection is performed to determine the maximum value of
% variation rate of the time-varying parameter.
%
% outputs: a table with the relaxation levels, values of gamma (maximum variation rate),
% number of scalar variables, number of LMI rows and computational times (seconds).
%
% date: 18/11/2008
% author: ricfow@dt.fee.unicamp.br
% 

pairs_g_d = [1 1;2 2;3 0;5 1;5 5;5 10;10 1;10 5;10 10];
mod_theta = 5; %|theta| <= 5
table = [];  
fprintf('running bissection ..\n');
for i=1:9
    g = pairs_g_d(i,1);
    d = pairs_g_d(i,2);
    gamma = bissection(mod_theta,g,d);
    [Ai,Bi,dot_alphas] = polytopic_representation(mod_theta,gamma);
    time = 0;
    for i=1:5        
        Info = synth_MOPB09_IJC_T3_c_yal(Ai,Bi,g,d,dot_alphas);
        time = time + Info.cpusec;
    end
    table = [table; [g d gamma Info.K Info.L time/5]];
end

fprintf('Final Results:\n');
fprintf('(g , d)\t gamma\t V\t  L\t    Time\n');
fprintf('-----------------------------\n');
for i=1:9
    fprintf('(%d,%d)\t %.3f\t %d\t  %d\t %.2f\n',table(i,1),table(i,2),table(i,3),table(i,4),table(i,5),table(i,6));
end

%_________________________________________________________________________________________
function gamma = bissection(mod_theta,g,d)


% bissection precision
maxgap = 0.001;

ub = 100;
lb = 0;

oncefeas = 0;
gap = ub - lb;
gamma = 0;
fprintf('[g=%d,d=%d] -> ',g,d);
while(gap > maxgap)

    rho = (ub + lb)/2;

    [Ai,Bi,dot_alphas] = polytopic_representation(mod_theta,rho);
    Info = synth_MOPB09_IJC_T3_c_yal(Ai,Bi,g,d,dot_alphas);

    if Info.feas == 1
        lb = rho;
        oncefeas = 1;
        gamma = rho;
        fprintf('%.3f ',rho);
    else
        if ~oncefeas
            lb = max(lb - 0.5,0);
        end
        ub = rho;
    end
    gap = ub-lb;    
end
fprintf('\n');


%_________________________________________________________________________________________
function [Ai,Bi,dot_alphas] = polytopic_representation(mod_theta,mod_dot_theta)

%affine representation
A0 = [-4.1 1;0 2];
A1 = [-3 0;-2 -3.2];
B0 = [3;2];
B1 = [0;-1];


%conversion to polytopic representation
Ai = [];
Bi = [];
theta=-mod_theta;
Ai = [Ai A0+theta*A1];
Bi = [Bi B0+theta*B1];

theta=mod_theta;
Ai = [Ai A0+theta*A1];
Bi = [Bi B0+theta*B1];

%conversion of the time-derivative bounds
t = mod_dot_theta/(2*mod_theta);
dot_alphas = [-t t;-t t];
