Problema

Este é o problema de controle ótimo mais próximo da realidade do dia-a-dia pois se propõe a solucionar um caso de otimização onde controles e estados possuem restrições. O que se deseja é minimizar a energia gasta do sistema sujeito à uma série de restrições. A técnica de controle ótimo para se resolver esse problema se chama técnica de Dreyfus, onde são geradas condições de entrada e saída nas restrições. Essas condições consistem na derivada da restrição (q) vezes até o aparecimento da variável de controle.
                Em termos matemáticos, o que se deseja é


sujeito à


vínculo dinâmico: 



onde o (.) representa a derivada temporal d(.)/dt.


vínculo de contorno:



vínculo do estado:


  (primeiro caso)
     (segundo  caso)


Solução


                Utiliza-se o funcional na forma de Mayer que consiste em transformar a integral em derivada e incluí-lo como variável de estado


    com


                Logo, o problema se torna um problema de valor de contorno em dois pontos, cujo sistema dinâmico completo será:


com as condições de contorno



Hamiltoniana



A expressão S(q) é encontrada pela derivada de q-vezes da restrição de estado até o aparecimento da variável de controle ( técnica de Dreyfus). No caso em questão,



Logo o valor de q é dois pois são necessárias duas derivadas em cima da restrição de estado para o aparecimento da variável de estado. Voltando a expressão da hamiltoniana



Equações Adjuntas (co-estado)



Mas pela condição final para a adjunta



Então, como a derivada da adjunta 3 é nula e seu valor final é 1, conclui-se que a adjunta do estado 3 (índice de performance) é constante e poderá ser substituída por 1 na equação da Hamiltoniana. Logo,



 Lei de controle ótimo


Ou seja,



A atuação do controle deverá ser analisada sob dois aspectos. Primeiro o caso onde o estado não está dentro das restrições, ou seja, o estado está livre. A segunda condição do controle ótimo deverá ser analisada quando o estado viola a restrição e então o controle deverá atuar para mantê-lo dentro dos parâmetros estabelecidos no projeto.

 

 

INTERVALOS LIVRES


Nesse caso, a ponderação (ou peso) da restrição do estado μ não está ativo e então



O sistema dinâmico completo para o problema de valor de contorno fica



INTERVALOS COM RESTRIÇÃO
Para intervalos onde o estado ultrapassa as restrições o controle sofre um penalidade ponderada pelo peso μ. Então, primeiro sabe-se que para os pontos de contato com a restrição, a variação da equação do vínculo deve ser nula, para mantê-lo sobre a linha de contato, ou seja,



e ainda o controle deve ser ótimo, ou



Então as duas equações são utilizadas para encontrar o controle ótimo com o fator de peso de sua correção. Resolvendo as duas equações



encontra-se a lei de controle ótimo para o estado



e o fator



O sistema dinâmico nesse caso será


 

Programação - Solução Numérica


 A solução possível para os casos enunciados de restrição do estado são resolvidas numericamente. Nos atuais dias, um bom software para a resolução de problemas de valor de contorno em dois pontos é o bvp4c.m ("colocation method") do Matlab. A programação para seu uso pode ser encontrada abertamente na internet, inclusive no site da Mathworks com exemplos. Basta observar que neste caso de restrição, o comando IF ("se") deve ser utilizado dentro do código principal para separar os casos fora da restrição dos casos dentro da restrição de estado.


Exemplo:


                                               if

                                                Else
                      

                                               End


A solução obtida usando o bvp4c.m é apresentada nos gráficos a seguir. A solução com as linhas em vermelho são para o caso 1, onde o estado x1 dentro da sua dinâmica, é proibido de ficar abaixo do valor -0,2. A linha em preto é para o segundo caso onde x1 deve sempre ficar acima de -1. O primeiro gráfico é o plano de fase das duas variáveis de estado mostrando as trajetórias ótimas. Interessante observar o gráfico à direita do primeiro gráfico onde se vê  perfeitamente o "batente" de restrição sendo respeitado. A linha reta mostra o estado percorrendo em cima da restrição, indicando que o ótimo é ficar no limite das imposições de projeto. No último gráfico à direita se vê a atuação do controlador e seu salto quando o estado entra na linha do limite aceitável. Esse salto mantém o estado sob a restrição até a dinâmica voltar ao estado normal permitido pelo projeto.

A restrição de contole foi adotada para sua atuação ficar entre

 

 

O ALGORITMO DO PROBLEMA EM MATLAB

O problema foi resolvido seguindo o algoritmo abaixo, usado para o programa principal que chama a rotina bvp4c.m

 

%PROBLEMA DE minima energia com vinculo de estado


% MIN  INT(0.5 u^2)dt
%   suj dx1=x2
%       dx2=-x2+u
%      
%     x1(0) =0 (dado)  x2(0)=0 
%     x1(tf)=1         x2(tf)=10
%     tf fixo
%     com x1>= L


%============================================================================


clear all
global T
global x10
global x20
global L
 
 
 
H =[]; G =[]; x = []; J = [];
%===============================INICIALIZACAO==============================
T=2;
x10=0;
x20=0;
 
%===============================restriçao de estado ==============================
L=-0.2;
 
%============================================== ===========================
 
type dreyfusInit;
type dreyfusM;
type dreyfusRes;
 
%===========================================================================
 
solinit = bvpinit(linspace(0,1,101),@dreyfusInit);
options = bvpset('Stats','on','RelTol',1e-6);
 
sol = bvp4c(@dreyfusM2,@dreyfusRes,solinit,options);
 
t = sol.x;
x = sol.y;
 
 
for i=1:1:size(t,2)
%   
control=-x(5,i);
if x(1,i)>L
 
   mm1(i)=-x(5,i);
 
else
    mm1(i)=x(2,i);
   
end
   
 
end;     
 
t=t'; x=x';
tsol=t(1:length(t))*T; xsol=x(1:1:length(t),:); cont= mm1(1:1:length(t))';
 
sol=[tsol xsol cont];
 
save dados.txt sol -ascii
 
 
%++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ +++++++++
%++++++++++++++++ SEGUNDO CASO ++++++++++++++++++++++++++++++++++++++++++++


clear all
global T
global x10
global x20
global L
 
H =[]; G =[]; x = []; J = [];
%===============================INICIALIZACAO ================================
T=2;
x10=0;
x20=0;
 
%===============================restriçao de estado =================================
 
L=-1;
 
%============================================== =============================

 
type dreyfusInit;
type dreyfusM;
type dreyfusRes;
 
%=================================================== ========================
 
solinit = bvpinit(linspace(0,1,101),@dreyfusInit);
options = bvpset('Stats','on','RelTol',1e-6);
 
sol = bvp4c(@dreyfusM,@dreyfusRes,solinit,options);
 
t = sol.x;
x = sol.y;
 
 
for i=1:1:size(t,2)
%   
control=-x(5,i);
if x(1,i)>L
    if control<=-3
        mm1(i)=-3;
    end
    if control>-3 & control<=100
   mm1(i)=-x(5,i);
 
   end
   if control>100
 
       mm1(i)=100;
   end
else
    mm1(i)=x(2,i);
    if mm1(i)<=-3
        mm1(i)=-3;
    end
    if mm1(i)>100
        mm1(i)=100;
 
    end
end
   
 
end;     
 
 
load dados.txt
time=t*x(6,length(t));
 
subplot(221)
plot(x(1,:),x(2,:),'-k',dados(:,2),dados(:,3),'-r')
xlabel('x1(t)')
ylabel('x2(t)')
grid
subplot(222)
plot(time,x(1,:),'-k',dados(:,1),dados(:,2),'-r')
xlabel('tempo')
ylabel('x1(t)')
grid
subplot(223)
plot(time,x(2,:),'-k',dados(:,1),dados(:,3),'-r')
xlabel('tempo')
ylabel('x2(t)')
grid
subplot(224)
plot(time,mm1,'-k',dados(:,1),dados(:,8),'-r')
xlabel('tempo')
ylabel('controle u(t)')
grid

 

 

 

 

 

Voltar ao índice de otimização dinâmica

 

 

Caso7: Minimizar energia com retrição de estado e controle