function saida = estabRob_T14(A)
% function saida = estabRob_T14(A)
%
% Determina a estabilidade robusta de um sistema linear incerto usando o seguinte problema de
% factibilidade:
%
% [X1Ai+Ai'X1'   Pi-X1+Ai'X2']
%  Pi-X1'+X2*Ai     -X2-X2'  ] < 0,   Pi > 0,   i=1,...,N
%
% As LMIs são programadas usando o YALMIP  e podem ser resolvidas
% por qualquer resolvedor de LMIs compatível com o YALMIP (SeDuM é o padrão).
%
% entradas : A  -> vértices da matriz dinâmica do sistema
%
% saídas : saida.factivel -> 1 se estável, 0 caso contrário.
%          saida.cpusec   -> tempo de cpu para resolver as LMIs (segundos)
%          saida.cpusec_m -> tempo de cpu para programar as LMIs (segundos)
%          saida.V        -> número de variáveis escalares.
%          saida.L        -> number de linhas de LMIs.
%          saida.P        -> variaveis solução P_i (se existirem).
%
% Autor: ricfow@dt.fee.unicamp.br
% Data:  07/09/2010
%
% Exemplo:
% clear all;
% M0=[-2.4     -0.6     -1.7     3.1
%     0.7     -2.1     -2.6     -3.6
%     0.5     2.4     -5.0     -1.6
%     -0.6     2.9     -2.0     -0.6];
% M1=[1.1     -0.6     -0.3     -0.1
%     -0.8     0.2     -1.1     2.8
%     -1.9     0.8     -1.1     2.0
%     -2.4     -3.1     -3.7     -0.1];
% M2=[0.9     3.4     1.7     1.5
%     -3.4     -1.4     1.3     1.4
%     1.1     2.0     -1.5     -3.4
%     -0.4     0.5     2.3     1.5];
% M3=[-1.0     -1.4     -0.7     -0.7
%     2.1     0.6     -0.1     -2.1
%     0.4     -1.4     1.3     0.7
%     1.5     0.9     0.4     -0.5];
% rho=1.49;
% A{1}=M0+rho*M1;
% A{2}=M0+rho*M2;
% A{3}=M0+rho*M3;
% saida = estabRob_T14(A)


% Determina as dimensões do problema.
ordem    = size(A{1},1);
vertices = size(A,2);

% Tempo de programação das LMIs.
saida.cpusec_m = clock;

% Inicializa a variável que armazenará as LMIs do problema.
LMIs = set([]);

% Cria as variáveis do problema.
for i=1:vertices
    P{i} = sdpvar(ordem,ordem,'symmetric');
end
X1 = sdpvar(ordem,ordem,'full');
X2 = sdpvar(ordem,ordem,'full');

% Entra com as LMIs do problema.
for i=1:vertices
    T11 = X1*A{i}+A{i}'*X1';
    T12 = P{i}-X1+A{i}'*X2';
    T22 = -X2-X2';
    LMIs = LMIs + set(P{i} > 0) + set([T11 T12; T12' T22] < 0);
end

% Finaliza o tempo de programação.
saida.cpusec_m = etime(clock,saida.cpusec_m);

% Determina o número de linhas de LMIs
saida.L = 0;
for i=1:size(LMIs,1)
    saida.L = saida.L + size(LMIs{i},1);
end

% Determina o número de variáveis escalares.
saida.V = size(getvariables(LMIs),2);

% Finalmente, o problema de otimização é resolvido.
solucao = solvesdp(LMIs,[],sdpsettings('solver','sedumi'));
saida.cpusec = solucao.solvertime;

% Acessa o menor resíduo primal.
p=min(checkset(LMIs));
saida.delta = p;

saida.factivel = 0;
% Capturando a solução (se existir).
if p > 0
    saida.factivel = 1;
    for i=1:vertices
        saida.P{i}     = double(P{i});
    end
end