Solucoes+em+Matlab

Prévia do material em texto

Mecatron – Projetos e Consultoria Júnior
Solução de problemas
usando Matlab®
Daniel Augusto Pereira
Alysson Fernandes Mazoni
2006
1
Sumário
Sistemas lineares.................................................................................................................................. 3
Problema exemplo........................................................................................................................... 3
Espaço de estados............................................................................................................................ 4
Função de transferência................................................................................................................... 8
Conversão de formas de representação............................................................................................8
Pólos, zeros, freqüência natural e fator de amortecimento............................................................. 9
Reposta a entradas padronizadas................................................................................................... 10
Resposta em freqüência................................................................................................................. 13
A ferramenta LTIview................................................................................................................... 16
Controle linear básico.................................................................................................................... 17
Controle proporcional............................................................................................................... 18
Outras configurações de controle..............................................................................................22
Ferramenta SISOtool (SISO – Single input single output)............................................................22
Simulação de sistemas de controle no Simulink................................................................................ 25
Otimização..........................................................................................................................................34
Problema exemplo ........................................................................................................................ 35
Otimização multiobjetivo.............................................................................................................. 36
Problema exemplo.....................................................................................................................37
Computação simbólica....................................................................................................................... 39
Exemplo: sistemas lineares............................................................................................................40
Exemplo: equações diferenciais.....................................................................................................41
2
Sistemas lineares
Ao chamar um sistema físico de linear estamos impondo sobre ele a hipótese da validade do 
princípio da superposição, ou seja, para uma combinação linear das entradas, a saída é uma 
combinação linear das saídas para cada entrada individualmente. Matematicamente, um sistema é 
linear se suas equações de evolução dinâmica (ao longo do tempo) são difenrenciais lineares e 
ordinárias.
Um sem número de sistemas de importância na engenharia são tratados como lineares. A 
exemplo: circuitos analógicos, motores de corrente contínua, sistema de troca de calor e estruturas 
flexíveis (metálica); outros são considerados lineares por trechos, comportam-se como lineares em 
uma faixa restrita de operação, como estruturas não lineares, atuadores hidráulicos, entre outros.
Problema exemplo
Para utilizar as funções do Matlab apropriadas para simulação e controle de sistemas 
lineares, utilizaremos o exemplo bastante comum de controle de um motor de corrente contínua. A 
representação utilizada aqui é a usual de um conjunto de equações diferenciais com o tempo como 
variável dependente, o que configura o que se chama de sistema dinâmico em tempo contínuo.
Um esquema da interação entre as partes elétrica e mecânica do motor é mostrado na figura 
abaixo:
Para esse modelo são escritas as equações diferenciais:
J ˙c−K i a=0
L i˙aRaRsi aK =V
˙=
em que  é a posição angular do rotor,  , a velocidade; J é a inércia associada; K, a constante 
eletromecânica; L, a indutância da armadura, Ra , a resistência de armadura; Rs , resistência 
em série com a armadura; e ia é a corrente de armadura.
Os valores medidos para um motor em particular são mostrados na tabela:
3
Parâmetro Valor
J 3,4E-4 Kg.m2
L 5,03E-3 H
K 1,3E-2 N.m/A
c 6,7E-5 N.s/m
Ra 4,03 
Rs 1 
Espaço de estados
Como é discutido em cursos de introdução a sistemas lineares e modelagem, um sistema 
como esse admite a conversão para forma de espaço de estados, o que corresponde a escrevê-lo na 
forma de um conjunto de equações diferenciais de primeira ordem. Para um sistema linear, isso 
leva, em notação vetorial, a
x˙=A xB u
y=C xD u
em que x representa um vetor com as funções no tempo que são soluções do sistema de equações 
diferenciais e é chamado vetor de estado. y é o vetor com as funções do tempo que se deseja 
observar: as saídas. u, é o vetor com sinais impostos ao sistema descrito ao longo do tempo de 
observação, é o vetor de entradas. As matrizes A, B, C e D contêm a representação do sistema 
dinâmico.
Um resultado elementar da teoria de sistemas lineares é o de que existem infinitas matrizes 
dessa forma que representam a mesma relação entre entrada e saída (o mesmo sistema), bastando 
escolher vetores de estado diferentes, ajustados entre um e outro por transformações lineares.
No caso particular em que o vetor de estados é constituído de variáveis físicas, diz-se que se 
possui um modelo físico nas matrizes citadas.
Para o caso do motor de corrente contínua citado um modelo físico pode ser escrito como:
[ ˙i a˙˙ ]=[−RaR sL −KL 0KJ −cJ 0
0 1 0
][ ia ][ 1L00 ]V
y=[0 1 0 ][ ia ][0 ]V
Dessa forma tem-se a entrada representada por V, a saída pela velocidade angular e o vetor 
de estados dados pela seqüência corrente, velocidade e posição. A saída observada, a velocidade 
angular, é definida pela forma como foi escrita a matriz C.
No Matlab, sugere-se declarar o modelo como:
4
>> J = 3.4e-4;
>> c = 6.7e-5;
>> K = 1.3e-2;
>> L = 5.03e-3;
>> R = 5.03;
>> A = [-R/L -K/L 0; K/J -c/J 0; 0 1 0];
>> B = [1/L; 0; 0];
>> C = [0 1 0];
>> D =[0];
Como já dito, há outras formas possíveis de representar a mesma relação matemática entre a 
entrada e a saída do sistema, porém com interpretações diferentes para os estados. Qualquer uma 
dessas representações pode ser usada para calcular, a cada instante de tempo desejado, a evolução 
dos estados no tempo usando os método usuais de integração direta. Dadas as matrizes, isso seria 
conseguido definindo uma função e usando um método como o ode45. Para tanto, é necessário 
também definir uma entrada, por exemplo senoidal.
function u=entrada(t)
u=sin(t);
function z = f_motor(t,x);
global A B
u = entrada(t);
z = A*x+B*u;
>> tin = 0; % tempo inicial
>> tfin = 10; % tempo final
>> x0 = [0;0;0]; % condição inicial
>> global A B
>> [t,y] = ode45('f_motor',[tin tfin],x0);
Note que as matrizes, usadas recorrentemente dentro e fora da função, foram definidas como 
variáveis globais. A figura abaixo mostra as respostas temporais do sistema, ou seja, mostra a 
evolução no tempo da corrente de armadura, da velocidade angular e da posição angular. Observa-
se que foi calculada a evolução temporal de todos os estados. 
>> u = entrada(t);
>> figure, plot(t,y,t,u);>> legend('Corrente de armadura [A]' , 'Velocidade Angular [rad/s]' , 'Posição angular [rad]' , 
'Tensão de alimentação [V]');
5
O Matlab contém estruturas de dados específicas para lidar com sistemas lineares 
representados nesse padrão. Um sistema linear é declarado numa estrutura de Matlab que mantém a 
representação em espaço de estados com o comando:
>> p = ss(A,B,C,D)
em que p passa a ser uma variável representando um sistema linear.
A partir dessa estrutura de dados, existem comandos no Matlab específicos para o cálculo da 
resposta temporal de um sistema. No caso da estrutura p definida o seguinte código pode ser 
utlizado:
>> t = 0:0.1:10;
>> u = sin(t);
>> x0 = [0;0;0]; 
>> figure, lsim(p,u,t,x0)
Nota-se que é visualizada apenas a saída selecionada pela matriz C do modelo no espaço de 
estados. Caso seja de interesse calcular a evolução temporal de todos os estados, bastaria redefinir a 
matriz C, criar sua nova representação no espaço de estados e então utilizar novamente o comando 
lsim.
>> C2=[1 0 0;0 1 0;0 0 1]; % Cálculo de todas as possíveis saídas
>> p2 = ss(A,B,C2,D);
>> figure, lsim(p2,u,t,x0)
Verifica-se a nova resposta na figura com as três saídas. Assim, fica evidente a importância 
do conhecimento da posição de cada estado dentro do vetor de variáveis de estado para a correta 
observação da saída desejada. 
6
0 1 2 3 4 5 6 7 8 9 1 0
- 1 0
- 5
0
5
1 0
1 5
2 0
2 5
3 0
 
 
C o r r e n t e d e a r m a d u r a [ A ]
V e l o c i d a d e A n g u l a r [ r a d / s ]
P o s i ç ã o a n g u l a r [ r a d ]
T e n s ã o d e a l i m e n t a ç ã o [ V ]
7
0 1 2 3 4 5 6 7 8 9 1 0
- 8
- 6
- 4
- 2
0
2
4
6
8
1 0
1 2
L i n e a r S i m u l a t i o n R e s u l t s
T i m e ( s e c )
Am
pli
tu
de
- 1
- 0 . 5
0
0 . 5
1
To
: 
O
u
t(1
)
- 1 0
- 5
0
5
1 0
1 5
2 0
To
: 
O
u
t(2
)
0 1 2 3 4 5 6 7 8 9 1 0
- 2 0
- 1 0
0
1 0
2 0
3 0
4 0
To
: 
O
u
t(3
)
L i n e a r S i m u l a t i o n R e s u l t s
T i m e ( s e c )
Am
pl
itu
de
Os gráficos gerados a partir do comando lsim podem ser editados clicando-se o botão direito 
do mouse e selecionando as opções desejadas. 
Caso o interesse seja armazenar os valores no tempo da saída em um vetor, ao invés de 
traçar um gráfico, basta fazer o seguinte:
>> y = lsim(p,u,t,x0); 
Função de transferência
Outra maneira de representar um sistema linear é através de funções de transferência. A 
função de transferência é definida como a razão entre a transformada de Laplace da função de saída 
sobre a transformada de Laplace da função da entrada. Conforme estudado no cursos elementares 
de equações diferenciais, essa razão, para condições iniciais nulas, pode ser facilmente obtida das 
equações diferenciais do sistema dinâmico representado. Há uma função de transferência para cada 
par de entrada e saída do sistema.
Para o sistema do motor de corrente contínua a função de transferência entre a alimentação e 
a velocidade é escrita como
 s 
V  s
= K /JL
s² RaRsL  cJ s c RaR sJL  K
2
JL
.
A transformada de Laplace de um sistema linear e invariante no tempo é sempre uma razão 
de polinômios na variável de Laplace (s). Desse modo, funções de transferência são especificadas 
através dos coeficientes dos polinômios do numerador e do denominador. A função de transferência 
é uma maneira unívoca de representar o sistema. Em código de Matlab:
>> np = [K/(J*L)];
>> dp =[1 ((Ra+Rs)/L+c/J) c*(Ra+Rs)/(J*L)+K^2/(J*L)];
>> H = tf(np,dp)
Alternativamente, pode ser definida da seguinte forma:
>> s = tf('s');
>> H = (K/(J*L))/(s^2 + ((Ra+Rs)/L+c/J)*s + c*(Ra+Rs)/(J*L) + K^2/(J*L))
A variável criada (H) é também um sistema linear, porém representado na forma de função 
de transferência. Os comandos de manipulação, caracterização e simulação de sistemas lineares 
que serão vistos a seguir aplicam-se tanto a uma quanto a outra forma de representação.
Conversão de formas de representação
As formas de representação definidas pelos comandos ss e tf podem ser convertidas de uma 
para a outra:
>> [a,b,c,d] = tf2ss(np,dp);
>> sys_mat = ss(a,b,c,d);
8
>> sys_fisico = ss(A,B,C,D);
>> [num,den]=ss2tf(a,b,c,d);
>> H=tf(num,den)
>> [num,den]=ss2tf(A,B,C,D);
>> H=tf(num,den)
É importante observar que as matrizes a, b, c e d fornecidas pelo comando tf2ss não são 
iguais àquelas deduzidas a partir do modelo físico (A, B, C e D). Isso ocorre porque a 
respresentação de um sistema no espaço de estados não é única e a representação padrão do Matlab 
não é, e nem poderia ser, igual à representação baseada no modelo físico (intuitivamente mais 
óbvia). É possível checar que os sistemas descritos pelas difenentes matrizes no espaço de estados 
são os mesmos. Para isso, basta verificar que os autovalores da matriz A (e consequentemente os 
pólos do sistema) são os mesmos para ambas as representações:
>> eig(A)
>> eig(a)
Pólos, zeros, freqüência natural e fator de amortecimento
Os pólos e zeros de um sistema linear podem ser obtidos diretamente a partir da função de 
transferência do mesmo, calculando-se as raízes dos polinômios do denominador e do numerador, 
respectivamente:
>> polos = roots(dp)
>> zeros = roots(np)
Os pólos podem ser obtidos através da matriz A do sistema no espaço de estados pois são 
iguais aos autovalores dessa matriz: 
>> polos = eig(A)
Existe ainda um comando específico para verificar os pólos de um sistema:
>> pole(H)
>> pole(sys_mat)
Graficamente é possível verificar a posição dos pólos e zeros no plano complexo:
>> pzmap(H)
>> pzmap(sys_mat)
9
Conhecidos os pólos e os zeros de um dado sistema, além de seu ganho estático (ganho em 
regime permanente, ou seja, s = 0 na função de transferência), pode-se definir um outra forma de 
representação de sistemas lineares:
>> k=dcgain(H);
>> z=roots(np);
>> p=roots(dp);
>> zp = zpk (z,p,k)
Para um sistema definido em uma variável de uma certa estrutura (função de transferência 
ou espaço de estados) é possível calcular a freqüência natural e o fator de amortecimento associados 
a cada um dos pólos:
>> damp(H) 
>> damp(sys_fisico) 
>> [wn,z,polos]=damp(H) 
Reposta a entradas padronizadas
Uma vez com o sistema definido em uma variável, as repostas a entradas notáveis da teoria 
de sistemas lineares, usadas para caracterizar o compostamento do sistema, são obtidas 
imediatamente. Por exemplo, a resposta ao degrau:
>> step(p)
10
- 1 0 0 0 - 9 0 0 - 8 0 0 - 7 0 0 - 6 0 0 - 5 0 0 - 4 0 0 - 3 0 0 - 2 0 0 - 1 0 0 0
- 1
- 0 . 8
- 0 . 6
- 0 . 4
- 0 . 2
0
0 . 2
0 . 4
0 . 6
0 . 8
1
P o l e - Z e r o M a p
R e a l A x i s
Im
a
gin
a
ry
 
Ax
is
>> step(H)
Com esse comando é traçado o gráfico da resposta ao degrau do sistema. Uma das vantagens 
da utilização da função step é a possibilidade de avaliar diretamente no gráfico os parâmetros da 
resposta ao degrau (sobressinal, tempo de estabilização, tempo de subida e valor de regime). Para 
isso, basta pressionar o botão direito do mouse, selecionar a opção characteristics e em seguida 
selecionar os parâmetros desejados (Peak Response, Settling Time, Rise Rime, Steady State). 
Ou ainda, pode-se desejar obter os valores no tempo da resposta e não traçar o gráfico:
>> [y,t] = step(p);
>> [y,t] = step(H);
Assim é possível armazernar os valores desejados para cálculos futuros ou para traçar 
gráficos de forma mais conveniente. Outro exemplo clásssico é a resposta ao impulso:
>> impulse(p)
>> impulse(s)
>> [y,t] = impulse(p);
>> [y,t] = impulse(s);
O gráfico gerado pela função impulse também pode sereditado clicando-se o botão direito 
do mouse.
11
12
S t e p R e s p o n s e
T i m e ( s e c )
Am
pli
tu
de
0 2 4 6 8 1 0 1 2 1 4 1 6 1 8 2 0
0
5
1 0
1 5
2 0
2 5
3 0
S y s t e m : H
F i n a l V a l u e : 2 5 . 7
S y s t e m : H
S e t t l i n g T i m e ( s e c ) : 1 3 . 2
S y s t e m : H
R i s e T i m e ( s e c ) : 7 . 4 3
Im p u l s e R e s p o n s e
T i m e ( s e c )
Am
pli
tu
de
0 2 4 6 8 1 0 1 2 1 4 1 6 1 8 2 0
0
1
2
3
4
5
6
7
8
S y s t e m : H
P e a k a m p l i t u d e : 7 . 5 2
A t t i m e ( s e c ) : 0 . 0 3 7 3
S y s t e m : H
S e t t l i n g T i m e ( s e c ) : 1 3 . 3
Resposta em freqüência
O termo resposta em freqüência significa a resposta em regime permanente de um sistema a 
uma entrada senoidal. Nos métodos de resposta em freqüência varia-se a freqüência do sinal de 
entrada dentro de um certo intervalo e estuda-se a resposta resultante. Uma das vantagens dos 
métodos de resposta em freqüência é a possibilidade de se utilizar os dados obtidos diretamente a 
partir das medições em um certo sistema físico sem necessitar recorrer a modelos matemáticos. 
Os métodos de resposta em freqüência foram desenvolvidos no período entre 1930 e 1940 
por Nyquist, Bode, Nichols e muitos outros. Esses métodos são dos mais poderosos na teoria de 
controle convencional.
As respostas em regime permanente às entradas senoidais podem ser obtidas diretamente a 
partir da substituição de s por na jω função de transferência, onde ω é a freqüência. Dessa forma, a 
resposta em regime permanente é dada por:
H  j=Me j
onde M é o módulo da função de transferência senoidal (ou relação de amplitude entre a entrada e a 
saída) e φ é a defasagem (ou diferença de fase entre a entrada e a saída). Ressalta-se que qualquer 
número complexo pode ser escrito em função um módulo e uma fase.
O diagrama de Bode é constituído de dois gráficos: um é o gráfico do módulo em dB (20 log 
|H(jω)|) de uma função de transferência senoidal; o outro é o gráfico do ângulo de fase (em graus). 
Ambos são traçados em relação à freqüência em escala logarítmica. Ou seja, os gráficos são 
traçados a partir do cômputo dos valores do módulo e da fase para toda uma faixa de freqüências.
A principal vantagem de se utilizar o diagrama de Bode é que a multiplicação dos módulos 
pode ser convertida em soma. Além disso, as curvas podem ser traçadas a partir de aproximações 
assintóticas simples. A expansão da faixa de baixas freqüências pelo uso da escala logarítmica de 
freqüência é muito vantajosa, visto que as características dos sistemas em baixas freqüências, na 
prática, são as mais importantes.
É muito simples obter o diagrama de Bode no Matlab. O código em Matlab a seguir mostra 
duas formas de se obter graficamente o diagrama utilizando a função bode. 
>> figure, bode(H)
 
>> w = 0:0.01:10e4;
>> [mag,phase] = bode(H,w);
>> magdb = zeros(1,length(mag)); % ajuste de dimensão
>> phase_mod = zeros(1,length(phase)); % ajuste de dimensão
>> magdb(1,:) = 20*log10(mag(1,1,:));
>> phase_mod(1,:) = phase(1,1,:);
 
>> figure, subplot(2,1,1), semilogx(w,magdb), axis([0 10e4 -150 50])
>> subplot(2,1,2), semilogx(w,phase_mod)
13
14
- 1 5 0
- 1 0 0
- 5 0
0
5 0
M
ag
n
itu
de
 
(dB
)
1 0 - 2 1 0 - 1 1 0 0 1 0 1 1 0 2 1 0 3 1 0 4 1 0 5
- 1 8 0
- 1 3 5
- 9 0
- 4 5
0
Ph
as
e 
(de
g)
B o d e D i a g r a m
F r e q u e n c y ( r a d / s e c )
Para analisar a estabilidade de um sistema linear pode ser utilizado o diagrama de Nyquist, 
que é o lugar geométrico de H(jω) no plano complexo à medida que ω varia de -∞ a +∞. Νo 
gráfico é possível perceber a marcação do ponto -1 no eixo real. Esse ponto é fundamental para a 
análise da estabilidade do sistema. Sem entrar em detalhes, nesse caso percebe-se que o sistema é 
estável porque o ponto -1+j0 não é envolvido pela curva fechada que representa o lugar geométrico 
de H(jω). 
15
1 0 - 2 1 0 0 1 0 2 1 0 4
- 1 5 0
- 1 0 0
- 5 0
0
5 0
1 0 - 2 1 0 0 1 0 2 1 0 4 1 0 6
- 2 0 0
- 1 5 0
- 1 0 0
- 5 0
0
A ferramenta LTIview
As funcionalidades de análise de sistemas LIT (lineares invariantes no tempo) são 
centralizadas na ferramenta LTIview, que permite acessar diversos tipos de gráficos já vistos, como 
dagrama de Bode, de Nyquist, resposta ao degrau, diagrama de lugar das raízes e outras como carta 
de Nichols e diagrama de valor singular. Nada pode ser alterado no sistema, esse é apenas um 
gerenciador gráfico para as diversas abordagens de visualizar um sistema LIT.
É acessada por:
>> ltiview
Usando o menu Edit – Plot configurations escolhe-se quatro gráficos por janela e as opções 
de resposta seguintes na ordem: resposta ao degrau, resposta ao impulso, diagrama de bode em 
magnitude e diagrama de Nyquist. A configuração e as respostas obtidas são mostradas nas duas 
figuras abaixo.
16
- 5 0 5 1 0 1 5 2 0 2 5 3 0
- 1 5
- 1 0
- 5
0
5
1 0
1 5
N y q u i s t D i a g r a m
R e a l A x i s
Im
ag
in
a
ry
 
Ax
is
Controle linear básico
Um dos objetivos de se obter um modelo de um sistema linear é projetar um controlador ou 
compensador para que se obtenham algumas características particulares na sua resposta. O caso 
mais comum de objetivo de controle é o regulador: um sistema que tem o objetivo de levar a 
resposta do sistema a zero em um tempo finito diante de qualquer entrada, ou tão próximo disso 
quanto possível. Já o controlador chamado rastreador tem a função de fazer a saída do sistema 
acompanhar a sua entrada, da forma mais próxima possível, para qualquer entrada.
É comum denominar o sistema a ser controlado de planta e a configuração com planta e 
controlador de malha fechada.
A teoria de controle de sistemas lineares é extensa e ainda em desenvolvimento, porém, 
pode-se dizer que a configuração de controle mais comum, ou seja, a conexão lógica e física entre 
um sistema e o seu controlador mais recorrente, é a seguinte:
17
Essa configuração corresponde ao controle em realimentação de saída. O controlador é um 
sistema linear que dispõe da informação da saída da planta. A entrada da planta com o controlador 
passa a ser a diferença entre a entrada original e o sinal de controle.
Há uma infinidade de outras configurações possíveis, por exemplo:
em que o controlador recebe o sinal de erro, ou seja, a diferença entre a entrada e a saída da planta.
Controle proporcional
O primeiro controlador possível é o controlador proporcional, representado por uma 
constante multiplicativa na segunda configuração mostrada acima. Isto é, regula-se o controle como 
um múltiplo do sinal de erro.
Para ajustar o sinal o comportamento do sistema controlado pelos pólos que ele possui pode-
se usar o diagrama no lugar das raízes. Esse diagrama traça a localização dos pólos do sistema 
quando há um controlador proporcional. Variando-se o ganho desse controlador, os pólos variam; o 
gráfico resultante é o citado diagrama. Ou seja, o ganha proporcional altera a localização dos pólos 
no plano complexo.
Para o sistema do motor de corrente contínua que usamos como exemplo, o diagrama de 
lugar das raízes produzido pode ser usado para ajustar um ganho de 32,9 levando a pólos de malha 
fechada complexos e de parte real -500 e imaginária 15,8. Esse valor foi escolhido ao acaso apenas 
para ilustração, levando-se em conta obviamente que o ganho deve levar a pólos de malha fechada 
que mantenham o sistema estável (pólos no semiplano esquerdo do plano complexo).
>> rlocus(p)
18
Planta
Controlador
SaídaEntrada
PlantaControlador SaídaEntrada
Para se obter a malha fechada como variável usa-se o comando:
>> mf = feedback(p,32.9)Nesse caso, a variável recém criada é um sistema linear que representa o motor original com 
o controlador proporcional. Nota-se que esse controlador está inserido no ramo da realimentação 
(conforme a primeira figura dos diagramas de controle). Podem ser traçados sua resposta ao 
impulso e mapa de pólos e zeros para verificar as características do novo sistema com a inclusão do 
controlador projetado.
>> step(mf)
>> pzmap(mf)
19
A posição atingida ao final das variações é próxima de 0,03, distante do valor do degrau (1). 
Esse erro pode ser corrigido multiplicando a entrada do sistema por um valor de ajuste.
No mapa de pólos e zeros, pólos são representados por X e zeros por círculos. São 
visualizados para esse caso, três pólos (um na origem e um par complexo conjugado com parte real 
-500 e imaginária 16,9) e um zero na origem.
Verificou-se que há uma difenrença na parte imaginária, que foi definida anteriormente 
como sendo 15,8. Isso surgiu devido à forma do diagrama, que possui uma mudança abrupta nessa 
20
região, ou seja, é apenas uma erro devido à forma como é feita a interpolação dos pontos do gráfico.
Uma outra abordagem seria utilizar o modelo com a saída em posição e não em velocidade. 
O diagrama para esse caso, depois de aproximado por zoom:
Usando o ganho de 11,8 o sistema obtido possui o seguinte comportamento de resposta ao 
degrau.
21
A figura mostra um sistema de elevada freqüência natural (com oscilações rápidas) 
comparativamente com o primeiro caso.
Usando o clique direito pode-se marcar características padrão da resposta do sistema, como 
o tempo de estabilização ou o valor de regime permanente.
Outras configurações de controle
Outros tipos de controlador e configurações são usadas. O primeiro passo nesse sentido é 
supor um controlador que não é estático (um ganho) mas um sistema também dinâmico. Os mais 
comuns são os controladores PID e avanço-atraso. Seus médotos de projeto podem ser executados 
em Matlab usando os comandos citados e testando a resposta.
Ferramenta SISOtool (SISO – Single input single output)
O projeto de sistemas de controle pode fazer uso da janela de funções acessada por:
>> sisotool
Acessa-se o menu File – Import para obter a variável com os dados do sistema. Ou, já 
informando o sistema a ser analisado:
>> sisotool(H)
22
Na janela que se abre na parte superior direita há um diagrama de blocos indicando qual a 
estrutura de controle usada. Clicando sobre o botão FS as posições dos blocos são alteradas. Nota-se 
que são considerados os seguintes tipos de bloco: a planta (G), o sensor (H – caso precise ser levado 
em conta no modelo), um pré-filtro da entrada (F) e o controlador (C).
A função do sensor tem o objetivo de modelar a dinâmica do transdutor utilizado para medir 
a saída do sistema quando da modelagem de uma planta real. O pré-filtro simula a característica de 
limitação em banda de freqüência do sinal que se pode impor à planta e ao controlador. No caso 
aqui apresentado supõe-se que o sensor mede perfeitamente, sem atraso e sem ganho a saída e que a 
entrada não possuirá componentes em freqüência que não possam ser realmente impostas à planta e 
ao controlador. Ou seja, F = H = 1.
À esquerda mostra-se a forma do controlador a ser ajustada (current compensator). Abaixo, 
os três gráficos padrão são o diagrama do lugar das raízes e os dois gráficos do diagrama de Bode.
O SISOtool se aplica apenas a sistemas com uma entrada e uma saída (como diz o nome em 
inglês). Sua operação consiste em escolher uma forma no domínio de Laplace para o controlador e 
ajustar seus parâmetros variando características requeridas nos gráficos de resposta mostrados.
Usando o modelo do motor com saída em velocidade e a configuração de controle fornecida 
pela janela inicial, clica-se sobre a caixa de edição do controlador.
Alí pode-se escolher a função de transferência desejada acrescentado ou retirando pólos e 
zeros. Para adicionar pólos ou zeros reais ou imaginários basta clicar no botão correspondente e 
digitar o valor inicial desejado. Para excluir, marca-se a opção delete.
23
Sugere-se incluir um pólo e um zero reais e negativos. Essa configuração corresponde aos 
controladores em avanço ou em atraso. Ou seja, um controlador com a seguinte função de 
transferência:
C  s=K s z
s p .
Escolhida a estrutura do controlador (clicando em OK), a ferramenta permite alterar 
qualquer desses três parâmetros dinamicamente. Nos gráficos de resposta mostrados basta clicar e 
arrastar os elementos que se deseja alterar e observar ar mudanças correspondentes nos outros 
gráficos.
Pólos de malha fechada são representados por quadrados. Pólos do controlador por X e zeros 
do controlador por círculos. Todos esses elementos pode ser arrastados com o mouse e as alterações 
correspondentes serão efetuadas. A curva de resposta do diagrama de Bode à direita pode ser 
arrastada para cima ou para baixo, mudando a banda de passagem e o ganho do controlador.
Um exemplo de teste assim com o sistema exemplo pode levar ao gráfico mostrado abaixo:
24
Os pólos e zeros do controlador podem ser acrescntados ou retirados também usando o 
clique direito sobre um dos gráficos e add pole/zero ou delete pole/zero. Os casos de derivador 
(differentiator) ou integrador (integrator) presentes por exemplo em um controlador PID podem ser 
acrescentador diretamente.
Outros gráficos podem ser mostrados na janela do SISOtool usando o menu View e 
escolhendo o gráfico ou dado a ser mostrado. No menu Analysis pode-se acessar a resposta do 
sistema em malha fechada invocando o LTIview, por exemplo, a resposta ao degrau é obtida em 
uma nova janela em Analysis – Rejection of step disturbance.
Simulação de sistemas de controle no Simulink
Os sistemas de controle trabalhados anteriormente podem ser simulados através do 
Simulink. Nota-se que os controladores devem ser previamente projetados porque o Simulink é uma 
ferramenta de simulação e de suporte à análise.
O primeiro caso de simulação a ser abordado é o controle proporcional da velocidade 
angular do mesmo motor de corrente contínua que vem sendo utilizado. Nesse caso, há duas 
diferentes maneiras de se posicionar o controlador: no ramo direto (amplificando o sinal de erro que 
entra na planta) ou no ramo inverso (amplificando o sinal de saída enviado para realimentação). As 
25
duas hipóteses foram simuladas. Para tanto, utilizou-se o mesmo ganho proporcional (Kp) projetado 
anteriormente, ou seja, um ganho de 32,9. 
A figura a seguir mostra os dois esquemas de controle proporcional utilizados para 
simulação.
Para a simulação, é utilizada como entrada uma função do tipo degrau unitário. Esse tipo de 
entrada é gerada a partir do bloco Step disponível na opção Sources da biblioteca de blocos. Ao ser 
inserido na área de criação, o bloco pode ser editado clicando-se duas vezes sobre o mesmo. Aqui o 
valor do degrau é mantido unitário, mas caso fosse necessário poderia ser escolhido uma 
“amplitude” diferente de 1 e também poderia ser alterado o instante no qual o degrau é aplicado. 
26
O sistema sob controle (motor CC) é definido a partir de seu modelo de estados (matrizes A, 
B, C e D) utilizando-se o bloco State-Space encontrado na opção Continuous. É necessário duplo 
clique sobre a caixa inserida para que se defina as matrizes de estado desejadas. Pode-se digitar as 
matrizes diretamente nos espaços disponíveis ou apenas digitar nos espaços adequados as variáveis 
já calculadas e disponíveis na área de comando do Matlab, conforme mostra a figura a seguir.
O ganho proporcional Kp é representado pelo bloco Gain encontrado na opção Commonly 
Used Blocks.Para alterar o ganho basta clicar no bloco e fazer a alteração.
A visualização da resposta é feita através do bloco Scope retirado da opção Sinks. Clicando-
se no bloco aparece uma janela com um gráfico bidimensional. Salienta-se que é padrão do Matla 
mostrar apenas os últimos 5000 pontos calculados, o que pode ser alterado clicando-se no segundo 
ícone no canto superior esquerdo (parameters) da janela gráfica, conforme mostra a figura abaixo.
27
A simulação das duas configurações de controle proporcional, Kp no ramo direto e Kp no 
ramo inverso, são mostradas, respectivamente, nas duas figuras a seguir. Observa-se que os 
resultados são os mesmos obtidos nas simulações anteriores, quando não foi utilizado o Simulink.
28
Para realizar o controle de posição, o procedimento é o mesmo adotado para o controle de 
velocidade. Para tanto, basta alterar a matriz C do sistema no espaço de estados para [0 0 1], o que 
pode ser feito diretamente na área de trabalho ou no próprio bloco State-Space. Utiliza-se o mesmo 
ganho proporcional projetado anteriormente, ou seja Kp = 11,8. Novamente são utilizadas as 
configurações com Kp no ramo direto e Kp no ramo inverso. Os diagramas são mostrados nas 
figuras a seguir.
29
A simulação dos dois sistemas de controle pode ser vista nas duas figuras a seguir. A 
primeira figura é resultado da aplicação de Kp no ramo direto e a segunda da aplicação de Kp no 
ramo inverso.
Outro tipo de controlador projetado anteriormente, com o auxílio da ferramenta SISOtool, 
foi um controlador em avanço (ou em atraso) de fase. Para inserir um controlador desse tipo, basta 
escolher a estrutura Zero-Pole encontrada na opção Continuous, e então definir os pólos, os zeros e 
o ganho desejados para realizar o controle em avanço (ou atraso). A figura a seguir mostra o 
esquema desse tipo de controlador.
A figura a seguir mostra a resposta ao degrau visualizada no bloco Scope.
Uma consideração importante a se fazer: é possível alterar o método (solver) utilizado pelo 
Simulink para calcular as respostas temporais que são mostradas nos blocos Scope. Por exemplo, 
pode-se mudar de passo fixo para passo variável, o que evita problemas quando há singularidades. 
Para isso, é necessário clicar na opção Simulation na barra de ferramentas superior e escolher a 
opção Configuration Parameters... e então aparecerá a janela mostrada a seguir, onde é possível 
fazer todas as alterações desejadas.
30
Também é possível definir o sistema a ser simulado diretamente a partir de sua função de 
transferência. Para isso basta definir seu diagrama de blocos, o qual pode ser facilmente construído 
no Simulink, conforme mostrado na figura abaixo. Nesse caso não são mostrados os valores dos 
ganhos, os quais são tirados diretamente da função de transferência entre alimentação e velocidade 
(motor CC).
31
Nota-se novamente a existência de um bloco somador, Sum, o qual está disponível na opção 
Commonly Used Blocks. Esse bloco pode ser editado conforme os demais. Assim, é possível alterar 
o número de entradas e os sinais das mesmas.
Há ainda blocos integradores, que são os responsaveis por integrar a derivada da posição 
para fornecer a posição. Esses blocos estão disponíveis na opção Continuous com o nome 
Integrator.
É possível criar subsistemas selecionando certos blocos, clicando-se o botão direito e 
selecionando a opção Create Subsystem. Selecionando todo o subsistema correspondente ao motor 
de corrente contínua, será criado um subsistema conforme mostrado na figura a seguir.
32
Também é possível fazer um controle proporcional conforme os anteriores, bastando para 
isso utlizar o subsistema gerado a partir do diagrama de blocos, os ganhos e o somador para a 
realimentação. Um controlador proporcional (para a velocidade) inserido no ramo direto, com o 
mesmo Kp já projetado ( Kp = 32,9).
A figura com o resultado da simulação pode ser observada abaixo. Percebe-se que os 
33
resultados são os mesmos já obtidos quando o sistema foi definido com o bloco State-Space.
Otimização
O problema matemático da otimização é o de encontrar, para uma dada função, seu valor 
mínimo (ou máximo) e o ponto ou intervalo em que ele ocorre. Isso pode ser escrito como:
min
x
f  x ou maxx
f  x .
Determinar x para o qual f(x) é mínimo ou máximo. x é potencialmente uma variável 
vetorial, ou seja, f pode depender de várias variáveis. Para o caso de valor mínimo tem-se um 
problema de minimização; maximização para o caso contrário. As funções da Matlab e outros 
programas e pacotes em geral são escritas de modo a contemplar apenas um dos dois problemas 
(geralmente o de minimização). Isso porque eles são intercambiáveis: maximizar f(x) é o mesmo 
que minimizar -f(x). A função f é chamada função objetivo.
Essa descrição é a de um problema de otimização irrestrita, ou seja, a variável x pode 
assumir qualquer valor. O conjunto de problemas mais complexo é o da otimização restrita, dado 
por:
min
x
f  x sujeito a g x 0h x=0
.
O que quer dizer que a variável (chamada em teoria de otimização de variável de decisão) 
deve pertencer a um conjunto delimitado por funções na forma de igualdades e desigualdades (f e g 
são funções vetoriais, podendo haver várias igualdades e várias desigualdades, que podem ser 
estritas ou não). Para o caso de desigualdades no sentido positivo, basta mutilplicá-las por -1, pois 
g(x) < 0 é o mesmo que -g(x) > 0.
Esses dois problemas representam a primeira divisão na otimização, há métodos específicos 
para problemas irrestritos e para problemas restritos (muito mais complicados).
Outros casos que representam divisores de águas na otimização são:
34
• função objetivo linear e restrições lineares – problema de programação linear, resolvido 
usando o algoritmo simplex ou um método de ponto interior (para muitas variáveis, 
centenas);
• função objetivo quadrática e restrições lineares – algoritmo de programação quadrática;
• função objetivo não linear e restrições não lineares – programação quadrática seqüêncial, 
método da penalidade, lagrangiano aumentado, métodos de ponto interior;
• função objetivo não linear, sem restrições – método de Newton;
• função objetivo descontínua ou “caixa preta” (sem informação das derivadas nem 
aproximações) – método de Nelder-Mead (busca), algoritmo genéticos;
• variável de decisão restrita a ser um número inteiro – programação inteira; no caso de a 
variável além de inteira ser binária tem-se a programação inteira binária.
O Matlab possui implementações para todas essas classes de problemas, em alguns casos, as 
funções identificam a ordem de grandeza do problema e trocam de algoritmo ou modo de operação 
dependendo do tamanho do problema.
Um problema historicamente anterior ao de otimização porém que pode ser tratado por 
algoritmos semelhantes é o de resolver um sistema de equações não lineares, ou seja, encontrar x tal 
que h(x) = 0. Por essa razão são incluídos no mesmo pacote de otimização no Matlab.
A funções apropriadas para lidar com esses problemas são mostradas na tabela abaixo. 
Todavia, não esgotam o universo de funções de otimização do Matlab:
Classe de problema Função(ões)
Minimização escalar (uma função e uma 
variável)
fminbnd
Minimização irrestrita fminunc, 
fminsearch
Programação linear linprog
Progração quadrática quadprog
Minimização restrita fmincon
Programação inteira binária bintprog
Solução de equação não linear fzero
Solução de sistema não linear fsolve
As funções de otimização exigem que sejam passadas a função objetivo e as funções da 
restrições. Isso pode ser feito de duas formas, definindo cada função na forma de um arquivo de 
Matlab (.m),definindo um função inline ou uma função anônima.
Problema exemplo 
Como exemplo propõe-se o seguinte problema de otimização restrita:
min f  x=ex14x1
22x2
24x1 x22x21
 sujeito a 
x1 x2−x1− x21,50
x1 x2−10
.
Convertendo a segunda restrição para a forma menor do que zero, podem-se definir duas 
funções no Matlab para representar a função objetivo e as restrições.
35
A função das restrições, pelo padrão da método a ser usado (fmincon) deve fornecer 
também restrições de igualdade como matriz vazia caso não existam.
function [c, ceq] = confun(x)
% Restrições de desigualdade.
c = [1.5 + x(1)*x(2) - x(1) - x(2); 
 -x(1)*x(2) - 10];
% Restrições de igualdade.
Ceq = [];
function z = fobj(x);
% Função objetivo.
Z = exp(x(1))*(4*x(1)^2+2*x(2)^2+4*x(1)*x(2)+2*x(2)+1);
Esse é um problema de duas variáveis. Uma chamada possível para a função fmincon, 
supondo um ponto inicial de busca da solução [-1, 1] é
x0 = [-1,1];
options = optimset('LargeScale','off');
[x, fval] = ... 
fmincon(@fobj,x0,[],[],[],[],[],[],@confun,options)
A variável options é usada para determinar parâmetros que controlam a execução do 
algorimo. Quase todas as funções do pacote de otimização admitem a estrutura com tais parâmetros. 
Essa estrutura pode ser criada pela função optimset informando nome e valor dos parâmetros 
que se deseja fazer diferentes dos padrões. Neste caso, escolheu-se usar a forma padrão do 
algoritmo e desativar o modo de problemas de grande porte, que é inicializidado por padrão e pode 
inclusive ser trocado automaticamente caso isso seja decidido pela função de otimização chamada.
O valor ótimo da função obtido é armazenado em fval e a variável de decisão ótima em x.
x = 
 -9.5474 1.0474 
fval =
 0.0236
A restrições foram quase exatamente satisfeitas (quase nulas), como pode ser verificado 
avaliando-as:
>> [c,ceq] = confun(x)
c =
 1.0e-007 *
 -0.9032
 0.9032
ceq =
 []
Otimização multiobjetivo
Todos esses problemas envolvem uma função objetivo escalar, ou seja, apenas um objetivo 
de otimização. Uma extensão da otimização está em imaginar que várias funções devem ser 
otimizadas simultaneamente. Nesse ponto entra a otimização multiobjetivo. A abordagem mais 
36
comum é a de definir uma função que combina de forma monotônica todas as funções envolvidas e 
minimizar essa função definida. Outra abordagem é a de obter computacionalmente o que se chama 
se superfície de Pareto: o conjunto dos pontos para os quais há pelo menos uma das funções 
objetivos otimizada sem que se possa melhorar outra função sem piorar alguma que já é ótimo.
O Matlab não implementa, até a presente versão, métodos de otimização que produzem 
automaticamente superfícies de Pareto, isso precisaria ser programado pelo usuário. No entanto, a 
primeira maneira citada de lidar com o problema pode ser usada.
Problema exemplo
O problema exemplo de minimização multiobjetivo consiste de encontrar valores ótimos 
para um controlador PID para otimizar o erro na resposta no controle de um sistema afetado por não 
linearidades de atuação.
Considere o sistema descrito pela função de transferência:
H  s= 1,5
50s343s23s1
,
em cuja entrada existam duas não linearidades: um limite de saturação e um limite no slew rate, que 
é o tempo mínimo para uma rampa de 1 V.
Coloca-se um controlador PID na entrada com parâmetros iniciais provisórios inicializados 
no espaço de trabalho e constrói-se o modelo no Simulink.
Coloca-se um degrau como entrada no sistema e o sinal de erro entre referência e saída é 
alimentado no controlador. Dentro do bloco do controlador colocam-se variáveis como valores dos 
parâmetros: Kp, Kd e Ki.
O diagrama em Simulink será então usado para simular repetidamente e obter a resposta do 
sistema variando os valores dos parâmtros do controlador até atingir valores ótimos no sentido de 
minimizar a soma dos quadrados dos erros entre a entrada em degrau e a saída. Aqui entra o 
conceito de multiobjetivo, pois a idéia é minimizar todos os erros simultaneamente, o critário 
adotador é o de minimizar a soma dos quadrados de todos eles.
Basta definir em Matlab uma função que invoca o diagrama e o simula dados os parâmetros 
37
do controlador. A sugestão aqui é a de usar uma função que produzirá os parâmetros ótimos pela 
seqüência: inicialização e chamada da função de minimização de quadrados.
Sugere-se também, por organização, definir a função que produz o sinal de erro a ser 
minimiado dentro da função geral.
function [Kp,Ki,Kd] = otimizaPID
% Otimização de parâmetros de um PID por mínimos quadrados.
 
pid0 = [0.63 0.0504 1.9688]; % valores iniciais
% Opções do algoritmo de otimização.
options = optimset('LargeScale','off','Display','iter',...
 'TolX',0.001,'TolFun',0.001);
% minimização do erro usando a funçõ definida abaixo
pid = lsqnonlin(@erro, pid0, [], [], options);
Kp = pid(1); Ki = pid(2); Kd = pid(3); 
 
 function F = erro(pid)
 % Produz o sinal de erro entre a saída e o degrau unitário.
 
 Kp = pid(1);
 Ki = pid(2);
 Kd = pid(3);
 % Parâmetros de simulação
 simopt = simset('solver','ode5','SrcWorkspace','Current');
 % Simulação
 [tout,xout,yout] = sim('optsim',[0 100],simopt);
 % Erro entre saída e degrau
 F = yout-1;
 
 end
end
O resultado obtido ao se invocar [Kp,Ki,Kd] = otimizaPID é
>> Kp
Kp =
 0.6300
>> Kd
Kd =
 1.9688
>> Ki
Ki =
 0.0504
38
Computação simbólica
A manipulação matemática de símbolos representando variáveis é chamada computação 
simbólica. É uma área da matemática ainda de desenvolvimento e cheia de problemas ainda não 
resolvidos, por exemplo, não há algoritmo definitivo ou mesmo robusto para simplificação de 
expressões matemáticas que opere para um conjunto amplo de tipos de expressões. Teoricamente, 
os programas de computação simbólica usam de estruturas de árvores para armazenam e manipular 
as expressões.
O Matlab implementa algumas técnicas de manipulação simbólica através do Symbolic 
Toolbox. É bastante conhecido que esse é o ponto fraco do programa e há pacotes mais adaptados 
para manipulação simbólica no mercado, como o Mathematica. No entanto a presença desse pacote 
atende a um requisito de completude do programa.
Define-se uma variável simbólica através de um nome de variável (como em outros casos) e 
uma string que é a representação escrita na tela desta variável. Para definir uma variável chamada 
var escrita na tela como x,
>> var = sym('x')
Para definir uma variável simbólica representada por uma string lexicamente idêntica ao seu 
nome pode-se usar o atalho (variável x escrita como x)
>> syms x
A utilizada da primeira versão do comando consiste em usar representações simbólicas 
familiares mesmo quando as letras usadas coincidem com nomes de outras variáveis em um 
programa longo.
As funções são baseadas em analogias com operações realizadas quando se realiza 
manipulação de expressões manualmente. Isto é, há funções para fatorar um termo, expandir 
produtos notáveis, entre outras na tabela abaixo:
Operação Função
Agrupar termos comuns collect
Expansão de polinômios expand
Fatorar factor
Separar numerador e 
denominador
numden
Procurar forma mais curta de 
escrever
simple
Simplificação simplify
Reescrever em termo de sub-
expressões
subexpr
Ainda, muitas funções numéricas são sobrecarregadas no pacote simbólico, pode ser 
aplicadas com a mesma chamada sobre variáveis numéricas ou simbólicas, outras são exclusivar do 
pacote, algumas mais comuns:
39
Operação Função
Derivada diff
Integral intMatriz jacobiana jacobian
Limite limit
Expressão para uma série symsum
Série de Taylor taylor
Determinante det
Autovalores e autovetores eig
Matriz inversa inv
Triângularização superior de 
matrizes
triu
Triângularização inferior de 
matrizes
tril
Valores singulares svd
Exponencial matricial expm
Composição de funções compose
Solução de equações solve
Função inversa finverse
Solução de equação diferencial dsolve
Gráfico de função simbólica ezplot
Ainda, todos os operadores do Matlab matém sua validade para variáveis simbólicas, tais 
como /, *, .*, ./, \,.
Exemplo: sistemas lineares
Definindo as variáveis simbólicas e o sistema:
>> syms a b c d r s
>> A = [a b; c d];
>> y = [r; s];
>> x = A\y
x =
 
 -(b*s-r*d)/(a*d-c*b)
 (a*s-c*r)/(a*d-c*b)
40
Exemplo: equações diferenciais
Para esse caso há um comando particular no pacote simbólico: dsolve. O operador de 
derivada é denotado por uma letra D seguida de um número representando a ordem da derivada. 
Para esta função, os argumentos são strings não havendo a necessidade de definir variáveis, porém, 
o resultado produzido é uma estrutura contento variáveis simbólicas.
>> resp = dsolve('2*Dy+3*y-x=5','Dx+x+(1/2)*y=0','y(0)=0','x(0)=0');
>> resp.x
 
ans =
 
-5/7+1/2*exp(-5/4*t)*(50/21*sin(1/4*3^(1/2)*t)*3^(1/2)+10/7*cos(1/4*3^(1/2)*t))
 
 
>> resp.y
 
ans =
 
10/7+exp(-5/4*t)*(20/21*sin(1/4*3^(1/2)*t)*3^(1/2)-10/7*cos(1/4*3^(1/2)*t))
Argumentos opcionais são as condições iniciais e a variável independente. Quando não 
informadas as respostas produzidas dependem de constantes padronizadas C1, C2, etc, e da variável 
t. Ou seja, outra maneira de chamar esta função seria:
>> resp = dsolve('2*Dy+3*y-x=5','Dx+x+(1/2)*y=0','x');
O que produziria respostas em função das constantes e com x como variável independente.
41
	Sistemas lineares
	Problema exemplo
	Espaço de estados
	Função de transferência
	Conversão de formas de representação
	Pólos, zeros, freqüência natural e fator de amortecimento
	Reposta a entradas padronizadas
	Resposta em freqüência
	A ferramenta LTIview
	Controle linear básico
	Controle proporcional
	Outras configurações de controle
	Ferramenta SISOtool (SISO – Single input single output)
	Simulação de sistemas de controle no Simulink
	Otimização
	Problema exemplo 
	Otimização multiobjetivo
	Problema exemplo
	Computação simbólica
	Exemplo: sistemas lineares
	Exemplo: equações diferenciais