Relatório Projeto Controlo de Processos G37 PDF

relatório PROCESS CONTROL Chemical Engineering MSc. Integrated Project for Final Evaluation Instituto Superior Técnico Professor Francisco Silva Lemos Professora Maria Amélia Lemos Maria Inês Nascimento – 113211 Hugo Marques – 113313...

relatório PROCESS CONTROL Chemical Engineering MSc. Integrated Project for Final Evaluation Instituto Superior Técnico Professor Francisco Silva Lemos Professora Maria Amélia Lemos Maria Inês Nascimento – 113211 Hugo Marques – 113313 Grupo 37 2025 Índice Resumo.................................................................................................................... 4 1. Introdução......................................................................................................... 5 2. Tanque de Pressão.............................................................................................. 6 2.1 Descrição do sistema.................................................................................... 6 2.2 Estudo do comportamento do sistema em Open Loop........................................ 7 2.3 Função de transferência para a variável manipulada............................................... 7 2.4 Função de transferência para a variável carga........................................................ 9 2.5 Perturbações sinusoidais para a variável manipulada e carga................................. 10 3. Tanque de nível................................................................................................ 12 3.1 Descrição do sistema....................................................................................... 12 3.2 Estudo do comportamento do sistema em Open Loop........................................... 13 3.3 Função de transferência para a variável manipulada............................................. 14 3.4 Função de transferência para a variável carga...................................................... 15 3.5 Perturbações sinusoidais para a variável manipulada e carga................................. 17 4. Estudo do comportamento dos sistemas em Closed Loop........................................ 19 4.1 Comportamento em Closed Loop para o Tanque de Pressão.................................. 19 4.1.1 Controlador PID....................................................................................... 20 4.1.2 Afinação do Controlador PID..................................................................... 23 4.1.3 Performance............................................................................................. 27 4.1.4 Estabilidade............................................................................................. 27 4.2 Comportamento em Closed Loop para o Tanque de nível...................................... 29 4.2.1 Controlador PID....................................................................................... 29 4.2.2 Afinação do controlador............................................................................ 31 4.3.3 Performance............................................................................................. 36 4.3.4 Estabilidade............................................................................................. 36 5. Aplicação de técnicas de reposta de frequências para análise de ambos os sistemas.... 38 5.2 Tanque de Pressão.......................................................................................... 39 5.3 Tanque de Nível............................................................................................. 49 6. Aplicação de modelos desenvolvidos para a análise e interpretação de resultados obtidos experimentalmente nas aulas de laboratório................................................................. 61 6.1 Tanque de pressão........................................................................................... 61 2  6.2 Tanque de nível.............................................................................................. 64 7. Conclusão....................................................................................................... 66 8. Referências...................................................................................................... 67 3 Resumo No presente projeto foi proposto desenvolver sistemas de controlo para dois sistemas distintos: um sistema de tanque de pressão e um sistema de tanque de nível. Para tal, este projeto dividiu-se em quatro etapas principais, cada uma com importância para o objetivo final de obter sistemas com um modelo de controlo e apropriado: 1. Estudo e Descrição dos sistemas em cadeia aberta Nesta primeira etapa, para cada um dos sistemas, foram descritas as equações de massa e determinadas as funções de transferência para as variáveis manipuladas e carga. Esta etapa teve como objetivo compreender a dinâmica básica dos sistemas sem a implementação de um controlador. 2. Estudo do sistema em cadeia fechada. Esta segunda etapa definiu-se pela implementação do sistema de controlo em cada um dos sistemas. O controlador usado para atuar sobre o sistema foi o controlador PID, constituído por diferentes modos operacionais. Realizou-se uma análise detalhada da estabilidade e do desempenho, utilizando diferentes métodos de controlo. Com isto a implementação do controlador recaiu sobre dois métodos: método de Ziegler- Nichols e o Controlo pelo Modelo Interno (IMC). Para analisar a estabilidade dos métodos, ambos foram estudados pela Matriz de Routh, que permite verificar a robustez dos controladores e pelo critério IAE, que mede a eficiência e a performance dos mesmos. 3. Análise da resposta do sistema à aplicação de frequências A terceira etapa consistiu num estudo das respostas dos sistemas quando submetidos a diferentes frequências, aprofundando a análise da estabilidade e desempenho dos controladores implementados. Para tal, foram utilizados gráficos de Bode e Nyquist, que permitem realizar uma avaliação mais abrangente das margens de estabilidade, da robustez e da resposta dos sistemas. 4. Ajuste com dados experimentais Por fim, os dados obtidos pela simulação do modelo teórico foram comparados com dados experimentais, de modo a validar o modelo e retirar conclusões sobre se o ajuste dos dados ao modelo teórico seria o mais apropriado. 4  1. Introdução Neste trabalho foram estudados dois tipos de sistemas: um Tanque de Pressão e um Tanque de Nível. O primeiro sistema analisado consiste num sistema de controlo concebido para a regulação de um tanque de pressão, um circuito contínuo de ar comprimido. (Fig. 1.) Figura 1 - Diagrama de Blocos para o Tanque de Pressão O segundo sistema analisado consiste num sistema de controlo para a regulação de um tanque de nível, um circuito de dois tanques ligados por uma tubagem, onde o caudal é controlado por uma válvula. (Fig. 2) O caudal entra no sistema, controlado por uma válvula de entrada, e sai do sistema, através de uma válvula de saída. Figura 2 - Diagrama de Blocos para o Tanque de nível 5 2. Tanque de Pressão 2.1 Descrição do sistema O sistema de controlo de pressão é constituído por um tanque de armazenamento onde o gás é mantido a uma pressão, p. Na entrada do tanque, existe uma válvula que controla o caudal de gás admitido, à pressão, pgrid, o ar comprimido é libertado com uma pressão (p), que, neste projeto, foi definida com um set-point (aplicado no controlador) de 1,9 atm. Na saída do tanque, encontra-se um sensor que está ligado ao início do sistema, influenciando assim a sua dinâmica. Existe também uma válvula que regula a quantidade de gás que sai do sistema, à pressão, patm. A variável manipulada é o caudal l1, que é ajustado através da abertura ou fecho das válvulas, permitindo assim estudar os seus efeitos no sistema, neste caso foi definida a 0,5 e para efeitos de carga foi o caudal de saída com a variação da abertura da válvula de saída, definida a 0,3. A variável regulada, que se pretende manter constante, é a pressão, p, de maneira a garantir que apenas a variável manipulada é alterada durante o estudo. As cargas presentes no sistema são a pressão de entrada, pgrid, a 2 atm, a pressão atmosférica patm a 1 atm, e os caudais l1 e l2. A pressão no interior do reservatório é controlada através de um mecanismo de retroalimentação que envolve um sensor, um controlador PID e um atuador. O sensor monitoriza continuamente a pressão de saída (p) e envia um sinal de retroalimentação ao controlador PID, que processa essa informação e gera um sinal de controlo. Este sinal é enviado ao atuador, que ajusta as condições de entrada de forma a garantir que a pressão no reservatório se mantém no nível desejado. Para este sistema, o balanço de massa pode ser descrito pela seguinte equação: 𝒅𝑷 𝑽 = 𝜶𝟏 𝒇(𝒍𝟏 ) 𝒑𝟐𝒈𝒓𝒊𝒅 − 𝒑𝟐 − 𝜶𝟐 𝒇(𝒍𝟐 ) 𝒑𝟐 − 𝒑𝟐𝒂𝒕𝒎 𝒅𝒕 𝑹𝑻 (1) Assumindo que as válvulas são lineares: 𝑭𝒆𝒆 𝟐 𝟐 𝒆𝒆 𝟐 𝟐 𝒆 = 𝜶𝟏 𝒍𝟏 𝒑𝒈𝒓𝒊𝒅 − 𝒑 = 𝑭𝒔 = 𝜶𝟐 𝒍𝟐 𝒑 − 𝒑𝒂𝒕𝒎 (2) 6  Utilizando as equações 1 e 2 é possível calcular o valor dos α (Tabela 1), assumindo que o caudal de ar é 120 L/min. Tabela 1 - Valores calculados para α. α1 α2 0,0064 0,0041 2.2 Estudo do comportamento do sistema em Open Loop Para estudar o comportamento do sistema em cadeia aberta recorreu-se ao método de Euler, um método numérico, em que é utilizado um passo de integração de dt = 0,1s e um atraso de θ = 1s, tornando assim possível a resolução das equações diferenciais. 𝒅𝑷 𝑷(𝒕 + ∆𝒕) = 𝑷(𝒕) + ∆𝒕 (3) 𝒅𝒕 Recorreu-se ainda a funções de transferência no espaço de LaPlace, para estudar as variações do sistema, seja na variável manipulada seja na carga. Assim, ao transformar equações diferenciais em equações algébricas mais simples, é, por isso, possível analisar a dinâmica do sistema, calculando a perturbação à saída através da soma das perturbações à entrada multiplicadas por uma função s. Como este sistema é de 1ª ordem, obtém-se a equação 4, como função de transferência de 1ª ordem com atraso, multiplicada pela perturbação. Onde K é o ganho do estado estacionário, o τ corresponde ao tempo que o sistema leva a reagir à perturbação, Δu corresponde à perturbação em degrau efetuada, e θ representa o atraso do sistema (delay). 𝑲∆𝒖 𝒚´(𝒔) = 𝒆𝜽𝒔 (4) (𝝉𝒔 + 𝟏)𝒔 Aplicando a expansão de Heaviside e a transformada inversa, obtém-se: 𝒕 𝜽 𝒚´(𝒔) = 𝑲∆𝒖 𝟏 − 𝒆 𝝉 ,𝒕 ≥ 𝜽 (5) 𝒚´(𝒕) = 𝟎, 𝒕 < 𝜽 (6) 𝒚(𝒕) = 𝒚´(𝒕) + 𝒚𝒆𝒆 (7) 2.3 Função de transferência para a variável manipulada Para obter as funções de transferência para a variável manipulada efetuou-se uma perturbação de 0,2. Através dos dados de simulação, arbitrou-se os parâmetros da função de 7 transferência, e, com recurso ao método do erro dos mínimos quadrados, utilizou-se a ferramenta Solver do Excel para minimizar o erro e obter os parâmetros reais. Figura 3 - Função de transferência para a variável manipulada com uma perturbação de 0,2. Tabela 2 - Parâmetros da função de transferência K (atm) τ (s) θ (s) 0,23 4,38 0,98 Obtidos agora os parâmetros, efetuaram-se as diferentes perturbações para avaliar a resposta do sistema às mesmas. Figura 4 - Função de transferência para diferentes perturbações na variável manipulada. Observando a figura 4, conseguimos concluir que, com os parâmetros apresentados na Tabela 2, a resposta do sistema aquando de uma perturbação de 0,2 é ideal. No entanto, quando as perturbações são negativas, os desvios entre a função de transferência e os dados utilizados 8 são maiores. Isto pode ser explicado pela calibração da função de transferência ter sido feita para perturbações positivas e por isso não conseguir ter uma resposta tão satisfatória. Além disso, estas funções, ao contrário do nosso sistema, são modelos lineares, o que explica o porquê das respostas para as perturbações positivas também apresentarem alguns desvios. Desta maneira, com a idealidade da resposta à perturbação, e com o sinal do ganho, podemos confirmar que os parâmetros se ajustam ao modelo utilizado. Visto que a perturbação é positiva, e que é uma ação direta, é de esperar que a pressão do sistema aumente, daí o valor positivo do ganho ser uma confirmação adicional do caráter positivo do modelo. 2.4 Função de transferência para a variável carga Tal como anteriormente, aplicou-se uma perturbação de 0,2 para a variável carga. Figura 5 - Função de transferência para a carga com uma perturbação de 0,2. Tabela 3 - Parâmetros da função de transferência para a carga. K (atm) τ (s) θ (s) -0,7 11 0,7 9  Figura 6 - Função de transferência para diferentes perturbações na carga. Neste caso, uma vez que a perturbação é feita na carga implica que quando a abertura da válvula de saída aumenta, a pressão no sistema vai diminuir e por isso o ganho do estado estacionário vai ser negativo. Assim, da mesma forma como a função de transferência da variável manipulada obtém desvios mais acentuados para as perturbações negativas, é previsível que para as perturbações positivas na carga, os desvios sejam, mais acentuados. Pelas mesmas razões que as anteriores, podemos confirmar que os parâmetros obtidos para a função de transferência ajustam-se ao modelo utilizado. 2.5 Perturbações sinusoidais para a variável manipulada e carga Após a aplicação de perturbações em degrau, foi realizada uma nova etapa experimental, utilizando perturbações sinusoidais, de maneira a caracterizar a dinâmica inicial do sistema. Este tipo de perturbação (Asin(ωt)) permite analisar a resposta em frequência, fornecendo informações cruciais sobre o comportamento do sistema em diferentes faixas de frequência, isto permite complementar a avaliação do ajuste dos parâmetros ao modelo. 𝑲𝑨 𝒕 𝜽 𝒚(𝒕) = 𝝎𝝉𝒆 𝝉 − 𝝎𝝉 𝐜𝐨𝐬 𝝎(𝒕 − 𝜽) + 𝐬𝐢𝐧 𝝎(𝒕 − 𝜽) (7) 𝝎𝟐 𝝉𝟐 + 𝟏 10  Utilizando os parâmetros das funções de transferência para a variável manipulada e para a carga, e aplicando várias perturbações, foram obtidos os gráficos representados nas figuras 7 e 8. Figura 7 - Resposta do sistema a perturbações sinusoidais (amplitude de 0,2), para a variável manipulada Figura 8 - Resposta do sistema a perturbações sinusoidais (amplitude de 0,2), para a carga Observou-se que, para ambas as variáveis, com diferentes amplitudes e frequências das perturbações, os desvios previstos pela função de transferência mantiveram-se próximos aos 11 valores experimentais, indicando que o modelo se ajusta bem às variações do sistema. Essa correspondência reforça a precisão dos parâmetros identificados previamente, confirmando a validade do modelo proposto para descrever a dinâmica do sistema. 3. Tanque de nível 3.1 Descrição do sistema O sistema de tanque de nível é um sistema composto por dois tanques interligados por uma tubagem, na qual o líquido que nela passa é controlada por uma válvula de permuta (p). Um dos tanques recebe o líquido que se pretende monitorizar pela abertura da válvula de entrada (e), a qual permite também variar e monitorizar diferentes quantidades de líquido. Para além das válvulas mencionadas, o sistema é constituído também por uma válvula de saída(s), que permite também o controle da quantidade de líquido dentro do tanque. Na análise deste sistema, o caudal de entrada foi considerado como a variável a manipular, de modo a avaliar a influência no comportamento na variável controlada, a altura do tanque 1, cujo valor inicial se estabeleceu como 20 cm. Para o resto da análise do sistema consideraram- se variáveis “carga” os caudais de saída e de permuta, visto que são variáveis que não sofrerão alterações “pensadas”, apenas por consequência de alterações na variável manipulada. Diferente ao observado para o sistema de tanque de pressão, o sistema de tanque de nível é definido por dois balanços de massa, referentes a cada um dos tanques presentes no sistema, constituídos pelas diferentes variáveis que neles atuam: 𝒅𝑽𝟏 = 𝑸𝒆(𝒍𝒆) − 𝑸𝒔(𝒍𝒔) − 𝑸𝒑(𝒍𝒑) (8) 𝒅𝒕 𝒅𝑽𝟐 = 𝑸𝒑(𝒍𝒑) (9) 𝒅𝒕 Para iniciar a análise do comportamento do sistema, estabeleceram-se alguns parâmetros:  As 3 válvulas (le, ls e lp) estão na mesma posição inicial (0,5)  Cada tanque teria uma área aproximada de 0,027 m2  Tanto o caudal de entrada, como o caudal de saída inicialmente teriam um valor de 0,0009 m3/min Para além destes parâmetros foi necessário realizar o cálculo de outros necessários à simulação do sistema: 12  Expressões Resultados 𝑸𝒆(𝒍𝒆) = 𝒍𝒆𝑸𝒎𝒂𝒙 𝑸𝒆𝒆 𝒆 𝑄𝑚𝑎x = 0,0018 m3/min 𝑸𝒎𝒂𝒙 = 𝒆𝒆 𝒍𝒆 𝑸𝒔 (𝒍𝒔 ) = 𝜶𝒔 (𝒍𝒔 )√𝒍𝒆 𝑸𝒆𝒆 𝒔 𝑸𝒆𝒆 𝒆 as = 0,004 𝒂𝒔 = 𝒆𝒆 = 𝒆𝒆 𝒍𝒆𝒆 𝒔 × 𝒉𝟏 𝒍𝒆𝒆 𝒆 × 𝒉𝟏 𝑸𝒔 (𝒍𝒔 ) = 𝜶𝒔 (𝒍𝒔 )√𝒍𝒆 𝑸𝒆𝒆 𝒔 𝑸𝒆𝒆 𝒆 ap = 0,004 𝒂𝒔 = = 𝒍𝒆𝒆 𝒔 × 𝒉𝒆𝒆 𝟏 𝒍𝒆𝒆 𝒆 × 𝒉𝒆𝒆 𝟏 Estes cálculos apenas são possíveis se considerarmos que as válvulas apresentam comportamento linear, ou seja que a sua abertura tem uma relação linear com o caudal de líquido que pode entrar no sistema. 3.2 Estudo do comportamento do sistema em Open Loop Para o estudo do sistema em Open Loop (Cadeia Aberta), também para o sistema de tanque de nível recorreu-se ao método de Euler que, como já referido anteriormente, se trata de um método numérico usado na resolução de equações diferenciais ordinárias, no qual é estabelecido um valor inicial de antemão. Para este sistema foi estabelecido um atraso de 5 minutos, associado a um tempo de integração de 0,2 minutos (20 segundos). Diferente ao observado para o sistema de tanque de pressão, para o sistema de tanque de nível foi obtida uma função de transferência de 2ª ordem, no estado de Laplace, com dinâmica de numerador (𝐾(𝑡 𝑠)) e com atraso associado (𝑒 ), visto que a manipulação é feita na variável de entrada e que se trata da associação de duas funções de transferência de 1ª ordem. 𝑲(𝒕𝒂 𝒔) ∆𝒖 𝜽𝒔 𝒚 𝒔(𝒕) = 𝒆 (10) (𝟏 + 𝒕𝟏 𝒔)(𝟏 + 𝒕𝟐 𝒔) 𝒔 De modo a interpretar o comportamento da variável no domínio do tempo, foi efetuada a transformada inversa, com recurso à expansão de Heaviside: 13  𝑡 −𝑡 ( ) 𝑡 −𝑡 ( ) (11) 𝑦 (𝑡) = 𝐾∆𝑢 1 + 𝑒 + 𝑒 ,𝑡 ≥ 𝜃 𝑡 −𝑡 𝑡 −𝑡 𝑦 (𝑡) = 0, 𝑡 < 𝜃 (12) 𝑦(𝑡) = 𝑦 (𝑡) + 𝑦 (13) 3.3 Função de transferência para a variável manipulada De modo a obter os parâmetros que caracterizam a função de transferência para o nosso sistema, foi efetuada uma perturbação em degrau (step disturbance) na variável manipulada de cerca de 0,1, significando que a válvula que controla o caudal de entrada abriu. Usualmente, a perturbação em degrau tem como intuito avaliar a robustez e a estabilidade do sistema a controlar. Com os dados recolhidos experimentalmente da altura com atraso, arbitrou-se parâmetros da função transferência e com recurso ao método do erro dos mínimos quadrados, minimizou- se este valor com a ferramenta Solver do Excel para se obter os parâmetros reais. 0,31 0,29 0,27 0,25 h1*(m) 0,23 0,21 0,19 0,17 0 50 100 150 200 t(min) h1atraso (m) h1atrasoFT* (m) Figura 9 - Função de transferência para a variável manipulada com uma perturbação de 0,1. Tabela 4 - Parâmetros da função de transferência para a variável manipulada. K (m) 0,88 Τa (min) 4,09 τ1 (min) 3,62 τ2 (min) 28,07 θ (min) 4,84 Os parâmetros obtidos para a função de transferência estão de acordo com o esperado, mais especificamente o sinal do ganho (K). Ao abrirmos mais a válvula de entrada, um maior 14 caudal de água passa pela válvula e mais água entra no tanque, resultando num aumento do nível de altura de líquido no tanque, daí o valor positivo do ganho. Como referido anteriormente, para além da perturbação feita para descobrir os parâmetros constituintes da função de transferência, foram feitas outras duas perturbações, para avaliar não só a estabilidade e robustez do sistema, mas para avaliar se os parâmetros da função demonstram um bom ajuste ao modelo experimental. 0,6 0,5 0,4 h1* (m) 0,3 0,2 0,1 0 0 50 100 150 200 time (min) h1*Δu=0,1 (m) h1*Δu=-0,1 (m) h1Δu=0,1FT* (m) h1*Δu=-0,1FT (m) h1*Δu=0,3 (m) h1Δu=0,3FT* (m) Figura 10 – Ajuste da função de transferência para diferentes perturbações na variável manipulada. Pela figura 10, é possível concluir que para uma perturbação (∆u) no le de 0,1, o ajuste é perfeito, o que é esperado visto que função de transferência foi obtida através desta mesma perturbação. Para as restantes perturbações, o ajuste também se enquadra, pois o desvio entre as curvas de cada uma das perturbações não é muito significativo. Outro detalhe que se evidencia da análise do gráfico é que a amplitude da perturbação pode ter influência sob o desvio entre a curva da função de transferência e a curva do atraso correspondentes à perturbação aplicada. 3.4 Função de transferência para a variável carga Equivalentemente ao que foi feito para variável manipulada, também para a válvula de saída (carga) foi aplicada uma perturbação, de modo a avaliar o comportamento da simulação. 15  0,25 0,2 0,15 h1*(m) 0,1 0,05 0 0 50 100 150 200 time (min) ∆Is=0,1 ∆Is=0,1 FT Figura 11- Função transferência para a carga com perturbação +0,1 Tabela 5 - Parâmetros da função de transferência para a carga K (m) -0,61 Τa (min) 10,43 τ1 (min) 12,22 τ2 (min) 15,86 Θ (min) 4,74 A abertura da válvula de saída promove uma maior saída de água do sistema, resultando numa diminuição da altura de líquido presente no tanque, logo o valor negativo de ganho. 0,35 0,3 0,25 0,2 h1*(m) 0,15 0,1 0,05 0 0 50 100 150 200 tempo (min) ∆Is=0,1 ∆Is=0,1 FT ∆Is=-0,1 ∆Is=-0,1 FT ∆Is=0,3 ∆Is=0,3 FT Figura 12 - Ajuste da função de transferência para diferentes perturbações na variável manipulada 16  De um modo geral a função de transferência aparenta ser um bom ajuste para o sistema, contudo são observáveis alguns desvios consideráveis para as restantes perturbações efetuadas. 3.5 Perturbações sinusoidais para a variável manipulada e carga Seguindo a mesma linha de análise feita para o tanque de pressão, após a perturbação em degrau foram realizadas diversas perturbações sinusoidais que, como já referido anteriormente, permitem estudar como o sistema responde a diferentes frequências. Para tal foi necessário associar à função de transferência que define sistemas de 2ª ordem parâmetros de amplitude e frequência: 𝐾(𝑡 𝑠) 𝐴𝑤 (14) 𝑦 𝑠(𝑡) = (1 + 𝑡 𝑠)(1 + 𝑡 𝑠) 𝑠 + 𝑤 Com recurso à expansão de Heaviside e à transformada inversa, foi possível obter a equação no espaço tempo: 𝑎1 ( ) 𝑎2 ( ) 𝑦 (𝑡) = 𝑒 + 𝑒 + 𝑎 sin 𝑤(𝑡 − 𝜃) + 𝑎 cos 𝑤(𝑡 − 𝜃) (15) 𝑡1 𝑡2 𝑡𝑎 𝐾𝐴𝑤 1 − 𝑡1 𝑎 = 1𝑡2 (16) 1 − 𝑡1 +𝑤 𝑡1 𝑡𝑎 𝐾𝐴𝑤 1 − 𝑡2 𝑎𝟐 = 𝑡1 1 (17) 1 − 𝑡2 +𝑤 𝑡2 (1 − 𝑡 𝑡 𝑤 ) + (𝑡 𝑡 )𝑡 𝑤 𝑎 = 𝐾𝐴 (18) (1 − 𝑡 𝑡 𝑤 ) + (𝑡 𝑡 ) 𝑤 ((1 − 𝑡 𝑡 𝑤 )𝑡 − (𝑡 𝑡 ))𝑤 𝑎 = 𝐾𝐴 (19) (1 − 𝑡 𝑡 𝑤 ) + (𝑡 𝑡 ) 𝑤 Tanto para a variável manipulada, como para a carga, foram realizados dois tipos de estudo: um onde se estabelecia um valor fixo de amplitude, fazendo variar a frequência; outro onde se efetuava a situação inversa. 17  0,27 0,25 0,23 h1*(m) 0,21 0,19 0,17 0,15 0 50 100 150 200 tempo (min) w=0,2 w=0,2-FT* (m) w=0,3 w=0,3FT* (m) w=0,5 w=0,5FT* (m) Figura 13 - Resposta do sistema a perturbações sinusoidais à variável manipulada, com amplitude 0,2. 0,29 0,27 0,25 0,23 h1*(m) 0,21 0,19 0,17 0,15 0 50 100 150 200 tempo (min) A=0,2 A=0,2FT* (m) A=0,1 A=0,1FT* (m) A=0,3 A=0,3FT* (m) Figura 14 - Resposta do sistema a perturbações sinusoidais à variável manipulada, com uma frequência a 0,2 De um modo geral, em ambos os estudos, os parâmetros da função de transferência parecem se ajustar com grande qualidade às perturbações sinusoidais efetuadas. Em adição, este tipo de perturbações permite também avaliar a robustez do sistema. Neste caso, apesar da perturbação, o sistema volta a atingir um estado estável, indicando um sistema robusto e equilibrado. 18  0,24 0,22 0,2 h1*(m) 0,18 0,16 0,14 0 50 100 150 200 tempo (s) ω=0,2 ω=0,2FT* (m) ω=0,3 ω=0,3FT* ω=0,5 ω=0,5FT* (m) Figura 15 - Resposta do sistema a perturbações sinusoidais à variável carga, com amplitude 0,2. 0,28 0,26 0,24 0,22 h1* (m) 0,2 0,18 0,16 0,14 0,12 0,1 0 50 100 150 200 tempo (min) A=0,2 A=0,2FT* (m) A=0,1 A=0,1FT* (m) A=0,3 A=0,3FT* (m) Figura 16 - Resposta do sistema a perturbações sinusoidais à variável manipulada, com uma frequência a 0,2 Nas perturbações realizadas à carga, a situação continua semelhante à encontrada para a variável manipulada: os parâmetros da função de transferência parecem se ajustar com grande qualidade às perturbações sinusoidais efetuadas. 4. Estudo do comportamento dos sistemas em Closed Loop 4.1 Comportamento em Closed Loop para o Tanque de Pressão Após a implementação do sistema em cadeia aberta, a segunda parte deste projeto recaiu sobre a implementação de um sistema de controlo. Ao implementarmos um controlador no sistema, o estudo do sistema passa a ser em cadeia fechada, onde, para este caso, vamos ter um 19 controlador que atua em feedback, ou seja, o controlador irá atuar consoante a perturbação sentida pelo sensor e consoante um ponto de ajuste (Set-Point) previamente definido. O objetivo principal do controlador é restaurar o sistema ao estado estacionário de forma segura, minimizando oscilações e garantindo estabilidade no processo Para este projeto, o controlador utilizado é um controlador PID, definido por 3 modos operacionais: Proporcional (P), Integral (I) e Derivativo (D). O modo proporcional, é responsável por reagir de forma imediata ao erro presente. Quanto maior o erro, maior será a correção aplicada pelo controlador, ou seja, quando o ganho do controlador aumenta, o erro do estado estacionário (offset) diminui. A ação integral acumula o erro ao longo do tempo e corrige a tendência de erro residual que a ação proporcional não conseguiu eliminar, diminuindo o tempo característico. A ação derivativa reage à taxa de variação do erro. Se o erro mudar rapidamente (como no caso de uma perturbação súbita), a ação derivativa antecipa essa mudança e aplica uma correção proporcional à velocidade de variação. O controlador é definido pela expressão 20, 𝒕 ∫ 𝜺𝒅𝒕 𝒅𝜺 𝑮𝒄 (𝒕) = 𝑲𝒄 𝜺+ 𝟎 + 𝝉𝑫 (20) 𝝉𝑰 𝒅𝒕 De modo a estudar a influência do controlador nos sistemas em estudo, foram feitas perturbações na carga e por tentativa erro chegou-se ao limiar da estabilidade, através dos parâmetros deduzidos. Posteriormente foram feitos ajustes ao controlador com recurso a dois tipos de controlador, o Ziegelr-Nichols e o Modelo Interno, tal como a análise de performance e de estabilidade. 4.1.1 Controlador PID Foi então utilizado um controlador PID para ajustar continuamente a entrada de pressão no sistema, ou seja, aplicado na variável manipulada, de maneira a minimizar o erro entre o valor desejado (set point) e o valor real medido da variável manipulada. A expressão 21 que descreve o controlador. 20  𝑴 = 𝑲𝒄 𝜺 = 𝑲𝒄 (𝑺𝑷 − 𝑪) (21) Apesar de não ser possível alcançar um controlador ideal devido à saturação dos sensores e atuadores e devido ao atraso dos equipamentos, a maneira a analisar a eficiência do controlador realizou-se com uma perturbação de 0,2 na carga (l2), tendo a abertura da válvula, no estado estacionário, permanecido nos 0,3. Posteriormente aumentou-se o ganho do controlador, Kc, até se alcançar o limiar da estabilidade, caracterizada por uma oscilação constante ao longo do tempo. Assim ao introduzir uma perturbação positiva na carga a pressão diminui no sistema, como a válvula de saída é aberta, e o objetivo do controlador é minimizar a influência dessa perturbação na pressão do sistema, a pressão vai tender a aproximar-se do valor de referência (set-point). Isto significa, que o controlador terá de corrigir o erro ao abrir a válvula de entrada (manipulada), no entanto devido ao aumento de pressão, a resposta será diminuir a abertura desta variável. Figura 17 – Determinação do limiar da estabilidade para a variação do Kc, com perturbação na carga de 0,2. Analisando o gráfico acima, podemos estabelecer o limiar da estabilidade como sendo Kc=35,4. Adicionando agora o modo integral, com o objetivo de eliminar o erro do estado estacionário (offset), a equação fica como a apresentada a seguir: 𝒕 𝟏 𝑴 = 𝑲𝒄 𝜺+ 𝜺𝒅𝒕 (22) 𝝉𝑰 𝟎 21  Assim quanto menor o τI, maior será a ação do modo integral. No entanto, as perturbações feitas no set-point podem levar a overshoots elevados e respostas lentas e instáveis, devido à saturação dos atuadores – Reset WIndup. De modo a tentar evitar este problema, é aplicado o anti-Reset Windup, com vista a suspender a ação integral enquanto o atuador está saturado. Para a parte integral, utilizou-se novamente o método de tentativa erro. Definiu-se o Kc, a utilizar como metade do valor obtido anterior, ou seja, 17,7. Figura 18 - Determinação do limiar da estabilidade para a variação do τI, com perturbação na carga de 0,2. Através da análise do gráfico anterior é possível estabelecer o limiar da estabilidade como τI = 1,05. No mesmo seguimento, a ação derivativa também apresenta as suas limitações – Derivative Kick. As perturbações no set-point podem levar a reações de forma exagerada. Isto porque, como a ação derivativa é proporcional à taxa de variação do erro, uma mudança repentina no set-point gera uma variação instantânea muito elevada no erro, o que resulta numa saída derivativa elevada e momentaneamente errática. Este comportamento pode causar oscilações e tornar o sistema instável. É então aplicado o anti-Derivative Kick de maneira a que a ação derivativa reaja apenas às variações da saída do sistema e não diretamente ao erro total, reduzindo as oscilações desnecessárias. 𝒕 𝟏 𝒅𝜺 𝑴 = 𝑲𝒄 𝜺 + 𝜺𝒅𝒕 + 𝝉𝑫 (23) 𝝉𝑰 𝟎 𝒅𝒕 22  𝒅𝜺 𝒅𝑺𝑷 𝒅𝑪 𝒅𝑺𝑷 = − ⟹ =𝟎 (24) 𝒅𝒕 𝒅𝒕 𝒅𝒕 𝒅𝒕 Este modelo deixa de ser aplicável se o sistema possuir muito ruído, visto que apenas o alongará. Para o estudo da parte integral, é necessário multiplicar o τI obtido no limiar da estabilidade por 3, como anteriormente, pelo método tentativa erro. Assim foi possível obter o gráfico abaixo representado, que demonstra um τD no limiar da estabilidade do sistema de 1,10. Figura 19 - Determinação do limiar da estabilidade para a variação do τD, com perturbação na carga de 0,2. Na tabela 6 estão então representados todos os parâmetros utilizados para o controlador. Tabela 6 - Parâmetros para o controlador PID Kc τI τD Limiar da estabilidade 35,4 1,05 1,10 𝐾𝑐 τ 3 τI Tentativa-Erro 2 3 17,17 3,15 0,37 4.1.2 Afinação do Controlador PID Apesar do método tentativa erro ser utilizado para modelar o controlador, levar o sistema até ao limite tantas vezes pode ser perigoso, tanto para o sistema como para os operadores do mesmo. Por esta razão, foi feita a afinação do controlador PID desenhado, pelos modelos de Ziegler-Nichols e Modelo Interno. 23  Ziegler-Nichols Utilizando o Kc deduzido anteriormente, agora Ku, K ultimate, foi medido o período entre oscilações, o resultado foi Pu=3,5. Tabela 7 – Parâmetros do controlador pelo método Ziegler-Nichols. Kc τI τD 𝑲𝒄 𝑃𝑢 𝑃𝑢 τ = τ = = 𝟎, 𝟔𝑲𝒖 2 2 10,62 1,75 0,44 Com estes parâmetros, foram feitas várias perturbações na carga e no set-point para posteriormente comparar com os resultados do modelo interno. Figura 20 – Resposta do modelo Ziegler-Nichols para perturbações no set-point. 24  Figura 21 - Resposta do modelo Ziegler-Nichols para perturbações na carga. Modelo Interno Para este modelo, foram utilizados os parâmetros da função de transferência da variável manipulada. Tabela 8 – Parâmetros do controlador pelo Modelo Interno Kc τI τD 𝜽 𝜏 𝟏 𝝉+𝟐 𝜃 τ = 𝑲𝒄 = τ = +𝜏 𝜃 𝑲𝝉 + 𝜽 2 22 +1 𝒄 𝟐 6,05 4,87 0,44 25  O valor de τc diz respeito a uma constante do modelo que Rivera et al. Recomenda que este seja simultaneamente maior que 0,8θ e 0,1τ. 2,05 2 Pressão (bar) 1,95 1,9 1,85 1,8 1,75 0 20 40 60 80 100 tempo (s) SP=0.1 SP=-0.1 SP=0.05 SP=-0.05 Figura 22 - Resposta do modelo interno para perturbações no set-point. 1,94 1,93 1,92 Pressão (bar) 1,91 1,9 1,89 1,88 1,87 1,86 0 20 40 60 80 100 tempo (s) l=0.2 l=0.05 l=-0.2 l=-0.05 Figura 23 - Resposta do modelo interno para perturbações na carga. Comparando ambos os métodos, é possível perceber que o Modelo Interno apresentou um desempenho superior ao método de Ziegler-Nichols. Este último apresenta uma resposta com um significativo overshoot, tanto para perturbações no set-point como para perturbações na carga. Embora o tempo de estabilização seja razoável, as oscilações iniciais indicam que a afinação por este método pode ser menos indicada para sistemas sensíveis ou que requerem alta estabilidade. 26  Contrariamente, o modelo interno demonstrou uma resposta mais suave, com menor overshoot e oscilações reduzidas. O tempo de estabilização também foi mais curto, especialmente para perturbações no set-point. Em suma, o Modelo Interno mostrou-se mais confiável, com melhor capacidade de rejeitar perturbações e manter a estabilidade do sistema. No entanto, o método Ziegler-Nichols pode ser útil como ponto de partida inicial para ajustes rápidos. 4.1.3 Performance Usualmente, podem ser realizados três tipos de critérios para avaliar a performance de um controlador, em que todos envolvem a realização do integral de erro. Esses três critérios são: o ITAE (integral do tempo ponderado ao erro absoluto), o ISE (Integral do erro ao quadrado) e o IAE (Integral absoluto do erro). Tendo em conta que o IAE é um critério intermédio e que envolve os outros critérios, foi o selecionado para avaliar a performance dos controladores face a perturbações na carga e no Set-Point. 𝐼𝐴𝐸 = |𝜀(𝑡)|𝑑𝑡 (25) Tabela 9 - IAE para o controlador de Ziegler-Nichols e Modelo Interno para várias perturbações no Set-Point e Carga SET POINT CARGA IAE SP=0,05 SP=-0,05 SP=-0,1 SP=0,2 SP=0,05 SP=-0,1 ZN 0,140607 0,172108 0,338365 0,081365 0,023048 0,054448 MI 0,269486 0,236859 0,474845 0,393873 0,100784 0,224594 Analisando os resultados do integral do erro absoluto, o controlador de ZieglerNichols mostra ter um menor erro em 4 das 6 perturbações efetuadas. 4.1.4 Estabilidade Relativamente à estabilidade dos controladores construídos foi tido em conta a Matriz de Routh, onde “um sistema é estável se e só se todos os polos da sua função transferência tiverem parte real negativa” e se um sistema sofreu alguma perturbação limitada no tempo, este apresenta também uma resposta limitada no tempo. 27  Para a matriz de Routh é condição necessária e suficiente que todos os zeros da equação característica tenham parte real negativa, e que os elementos mais à esquerda da matriz sejam positivos. Assim para o Closed loop temos a equação apresentada como denominador da seguinte função de transferência. 𝐺 𝐾 𝐺𝐺 𝐺 𝐶´ = 𝐿 + 𝑆𝑃′ (26) 1+𝐺 𝐺 𝐺 𝐺 1+𝐺 𝐺 𝐺 𝐺 𝐺𝑣 (função transferência do atuador) e 𝐺𝑚 (função transferência do sensor) são considerados com valor 1. 𝐺𝑐 (função transferência do controlador) e 𝐺𝑝 (função de transferência da manipulada) já foram anteriormente descritas. Simplificando o polinómio do denominador temos, 𝐾𝑝 1 𝑃(𝑠) = 1 + 𝑒 ×𝐾 1+ +𝜏 𝑠 (27) 1 + 𝜏𝑠 𝜏𝑠 Aplicando a aproximação de Padé, 𝜃 𝐾𝑝 1 − 2 𝑠 1 𝑃 (𝑠) = 1 + ×𝐾 1+ +𝜏 𝑠 (28) 1 + 𝜏𝑠 1 + 𝜃 𝑠 𝜏𝑠 2 Simplificando, 𝑎 𝑠 +𝑎 𝑠 +𝑎 𝑠+𝑎 =0 (29) Como, 1 1 𝑎 = − 𝐾 𝜏 𝜏 𝐾𝑝𝜃 + 𝜏 𝜏 𝜃 (30) 2 2 1 1 𝑎 = 𝐾 𝜏 𝜏 𝐾𝑝 − 𝐾 𝜏 𝐾𝑝𝜃 + 𝜏 𝜏 + 𝜏 𝜃 (31) 2 2 1 𝑎 = 𝐾 𝜏 𝐾𝑝 − 𝐾 𝐾𝑝𝜃 + 𝜏 (32) 2 𝑎 = 𝐾 𝐾𝑝 (33) 𝑎 𝑎 −𝑎 𝑎 𝑏 = (34) 𝑎 𝑏 𝑎 −𝑎 𝑎 𝑐 = =𝑎 (35) 𝑏 28  Obtemos a seguinte matriz de Routh 𝑎 𝑎 0 𝑎 𝑎 0 𝑏 0 0 𝑐 0 0 E fazemos os cálculos para ambos os controladores, Ziegler-Nichols 2,83 4,84 0 8,31 2,45 0 4,01 0 0 2,45 0 0 Modelo Interno 8,97 10,99 0 38,24 1,39 0 10,66 0 0 1,39 0 0 Assim como todos os valores da coluna da esquerda são positivos, ambos os controladores, obtidos pelos métodos de afinamento, são estáveis. 4.2 Comportamento em Closed Loop para o Tanque de nível Em semelhança ao realizado para o sistema de tanque de pressão, foi equacionado também para o sistema de tanque de nível a implementação de um sistema de controlo, que consoante eventuais perturbações efetuadas no sistema, conseguisse detetar e atuar da melhor forma, através de um sensor e atuador, respetivamente. A implementação do controlador é essencial para que, face a uma perturbação, o sistema esteja capacitado a retornar ao estado de estabilidade inicial. 4.2.1 Controlador PID Mais uma vez, o controlador PID foi o selecionado para, através dos seus 3 modos operacionais, atuar sob as válvulas das variáveis manipuladas. Tendo por referência a teoria e as equações presentes na tabela 6, através do método arcaico de tentativa e erro, foram calculados os 3 parâmetros referentes aos 3 modos operacionais do controlador. No caso do tanque de nível, seguindo o referido em , começou-se por efetuar uma perturbação na carga com valor de 0,2, o equivalente ao abrir a válvula de saída. Isto fez com que o estado estacionário alterasse o seu valor de 0,5 para 0,7. 29  Consequentemente, foram impostos os 3 modos operacionais do controlador, todos desligados, ou seja, com um valor de Kc=0, um valor de ti de 100000 e por fim um valor de tD=0. Postumamente, variando cada um dos parâmetros individualmente, obteve-se o seu limite de estabilidade face à perturbação aplicada. 0,22 0,2 0,18 h1* (m) 0,16 0,14 0,12 0,1 0,08 0 50 100 150 200 Tempo (min) kc=0 kc=5 kc=10 kc=11,05 Figura 24 - Determinação do limiar da estabilidade variando o Kc, com uma perturbação na carga de 0,2. 0,21 0,2 0,19 h1* (m) 0,18 0,17 0,16 0,15 0 50 100 150 200 Tempo (min) tau I =10 Tau I = 6 Tau I = 5,7 Tau I = 5,55 Figura 25 - Determinação do limiar da estabilidade variando o TI, com uma perturbação na carga de 0,2 30  0,205 0,2 0,195 0,19 h1* (m) 0,185 0,18 0,175 0,17 0,165 0,16 0 50 100 150 200 Tempo (min) Tau D = 0 Tau D =1 Tau D= 0,5 Tau D = 1,3 Figura 26 - Determinação do limiar da estabilidade variando o Tau D, com uma perturbação na carga de 0,2. Tabela 10 - Parâmetros do controlador PID obtidos pelo método de tentativa e erro Parâmetros do controlador PID Kc (m) 𝜏𝐼(𝑚𝑖𝑛) 𝜏𝐷(𝑚𝑖𝑛) Valor correspondente ao limiar 11,05 5,55 1,3 de estabilidade Valores obtidos para o método 𝐾 𝑡 = 5,525 3 × 𝑡 = 16,65 = 0,433 tentativa e erro 2 3 A análise das figuras 24, 25 e 26 permite afirmar que teoricamente o controlador atuou de forma correta, consoante o valor e a variável na qual a perturbação foi imposta: ao aplicarmos uma perturbação na válvula de saída, a altura de líquido no tanque – variável regulada- sofre um decréscimo. Logo o atuador tem de apresentar uma ação inversa para que o nível do tanque volte ao estado inicial e, para isso a válvula de entrada tem de abrir, de modo a que o caudal de entrada aumente também. Isto reflete-se a um ganho da variável de manipulada positivo, o qual deve ser também refletido no ganho do controlador (Kc), que é o que acontece. 4.2.2 Afinação do controlador Tal como referido anteriormente, é necessário assegurar que os parâmetros correspondentes aos modos operacionais do controlador estejam bem definidos para que o atuador possa, consequentemente, atuar da melhor forma. Apesar de se ter recorrido ao método de tentativa e erro, este não é o mais indicado. Posto isto, recorreu-se a dois métodos de afinação do controlador: Ziegler-Nicholas e o Modelo Interno. 31  Da mesma maneira que foi realizado para o sistema de tanque de pressão, para o sistema de tanque de nível a implementação do método de Ziegler-Nichols fez-se iniciar pela aplicação de uma perturbação no Set-Point de -0,05, seguida da variação do valor de Kc até se atingir o limiar de estabilidade. A partir da curva referente ao limiar de estabilidade, obteve-se o período oscilatório, pela distância entre dois picos, tornando-se possível o cálculo dos restantes parâmetros, através das equações presentes na Tabela 7. Tabela 11 - Parâmetros do controlador pelo método de Ziegler-Nichols, para perturbação no Set-Point de -0,05 𝑲𝒖 (𝒎𝒊𝒏) 𝑷𝒖 (𝒎𝒊𝒏) 𝑲𝒄 (𝒎𝒊𝒏) 𝒕𝑰 (𝒎𝒊𝒏) 𝒕𝑫 (𝒎𝒊𝒏) 10,6 16,2 6,36 8,1 2,03 Definidos os parâmetros para o controlador PID segundo o método de Ziegler-Nichols, foram realizadas perturbações adicionais no Set-Point e na carga, de modo a analisar o comportamento do sistema e a resposta do controlador. 0,4 0,35 0,3 0,25 h1*(m) 0,2 0,15 0,1 0,05 0 0 50 100 150 200 tempo (min) ∆SP =-0,05 ∆SP =-0,1 ∆Sp =0,1 ∆SP =0,05 Figura 27 - Resposta do controlador Ziegler-Nichols para diversas perturbações no Set-Point 0,24 0,23 0,22 0,21 h1*(m) 0,2 0,19 0,18 0,17 0,16 0 50 100 150 200 tempo (min) ∆LS =0,2 ∆LS =0,1 ∆LS =-0,1 ∆LS =-0,2 Figura 28 - Resposta do controlador Ziegler-Nichols para várias perturbações na carga 32  Consequentemente, por motivos de comparação, foi implementado o método do modelo interno, o qual se baseia nos parâmetros já determinados para a função transferência da variável manipulada. 1𝑡 𝑡 −𝑡 𝐾𝑐 = (36) 𝐾 𝑡 +𝜃 𝑡 =𝑡 +𝑡 −𝑡 (37) 𝑡 𝑡 − 𝑡 (𝑡 + 𝑡 − 𝑡 ) 𝑡 = (38) 𝑡 +𝑡 −𝑡 Tendo por base as equações 36, 37 e 38, referentes a um sistema de segunda ordem com atraso e dinâmica de numerador, calcularam-se os parâmetros do controlador. Contudo, para 𝑡 obteve-se um valor negativo, significando que este modelo não seria o mais apropriado para a função de transferência obtida. Como tal, assumiu-se que, devido ao facto de que função de transferência para variável manipulada continha valores de 𝑡 e 𝑡 muito similares, seria possível ajustar a função de transferência de segunda ordem, para uma função de 1ªordem com atraso. Tabela 12 - Parâmetros da função transferência de 2º Ordem e da nova função transferência de 1º Ordem. Função de Transferência de 2ª ordem Função de Transferência de 1ª ordem K (m) 0,88 K (m) 0,88 𝒕𝒂 (𝐦𝐢𝐧) 4,09 𝑡(min) 27,82 𝒕𝟏 (𝐦𝐢𝐧) 3,62 𝜃(min) 4,54 𝒕𝟐 (𝐦𝐢𝐧) 28,07 𝜽 (𝐦𝐢𝐧) 4,84 Por comparação dos valores obtidos para ao parâmetro em ambas as funções de transferência, é possível concluir que não existe uma discrepância significativa entre parâmetros, pelo que pode ser considerado, que a aproximação é um ajuste de boa qualidade. 33  0,3 0,27 0,24 h1*(m) 0,21 0,18 0,15 0 50 100 150 200 tempo (min) h1* h1* FT Figura 29 - Ajuste da nova função transferência de 1º Ordem, obtida para uma perturbação de 0,1 na variável manipulada. A figura 29 demonstra que, mesmo com a aproximação da função de transferência para uma equação de primeira ordem, o sistema, após a perturbação, consegue chegar a um novo estado estacionário. Posto isto, foi possível a obtenção dos parâmetros do controlador segundo o método do modelo interno, segundo as equações 39, 40 e 41. 𝜃 1 𝑡+2 𝐾𝑐 = 𝐾𝑡 +𝜃 (39) 2 𝜃 𝑡 = +𝑡 (40) 2 𝑡 𝑡 = 𝑡 2𝜃 + 1 (41) O valor de 𝑡 foi escolhido segundo o referido por Rivera et al. – 𝑡 deve ser superior a 0,8𝜃 e a 0,1 𝑡. Como tal, estipulou-se que 𝑡 = 5 e em consequência, obtiveram -se os parâmetros para o controlador PID. Tabela 13 - Parâmetros do controlador PID pelo método do Modelo Interno 𝑲𝒄 (𝒎) 4,71 𝒕𝑰 (𝐦𝐢𝐧) 30,09 𝒕𝑫 (𝐦𝐢𝐧) 2,10 𝒕𝒄 (𝐦𝐢𝐧) 5 34  Definidos os parâmetros para o controlador PID segundo o método do Modelo Interno, foram realizadas perturbações adicionais no Set-Point e na carga, de modo a analisar o comportamento do sistema e a resposta do controlador. 0,35 0,3 0,25 h1*(m) 0,2 0,15 0,1 0,05 0 0 50 100 150 200 tempo(min) ∆SP=0,05 ∆SP=0,1 ∆SP=-0,1 ∆SP=-0,05 Figura 30 - Resposta do controlador Modelo Interno para várias perturbações no Set-Point 0,3 0,25 h1*(m) 0,2 0,15 0,1 0 50 100 150 200 tempo(min) ∆Is=0,1 ∆Is=0,3 ∆Is=-0,1 ∆Is=-0,3 Figura 31 - Resposta do controlador Modelo Interno para várias perturbações na carga Analisando o comportamento dos dois métodos perante as duas perturbações, é possível concluir que:  O controlador obtido segundo o método de Ziegler-Nichols parece conseguir retornar ao estado estacionário mais rapidamente do que o controlador obtido segundo o método do modelo interno, aquando de perturbações na carga.  Para perturbações no Set-point, a situação é inversa. O controlador segundo o modelo interno apresenta respostas mais rápidas e com menos presença de ruído e oscilação. 35  Tendo isto em consideração, é difícil afirmar que um método seja mais apropriado que outro, pelo que é necessário ter mais fatores em consideração como a estabilidade ou performance. 4.3.3 Performance Existem diversos métodos designados a avaliar a performance de um controlador. Para manter a concordância e a similaridade no trabalho, mantivemos a realização da análise da performance pelo método IAE, o qual se caracteriza por ser um método intermédio dos métodos ITAE e ISE. Assim sendo o método IAE consegue avaliar tendo por base mais tipos de erros. Para tal foram realizadas perturbações no Set-Point e na carga, de modo a perceber qual controlador apresentava menor quantidade de erros. Tabela 14 - IAE para o controlador de Ziegler-Nichols e Modelo Interno para várias perturbações no Set-Point e Carga. Set Point Carga IAE ∆SP=0,05 ∆SP=-0,05 ∆SP=-0,1 ∆Is=0,05 ∆Is=-0,1 ∆Is=0,2 ZN 1,08255 1,07792 1,81175 0,60881 0,64553 0,69219 MI 0,56920 0,56331 1,02491 0,39560 0,70549 1,31683 Conclui-se que o controlador obtido segundo o método do modelo interno parece ter mais performance global, visto que contém um menor erro em quatro das seis perturbações. 4.3.4 Estabilidade Para avaliar a estabilidade dos métodos designados para a construção do nosso controlador, mais uma vez recorreu-se à matriz de Routh, a qual assenta na premissa que: para um sistema ser considerado estável, necessariamente para todos os zeros da equação característica terem parte real negativa, os elementos mais à esquerda da matriz serem positivos. Como referido anteriormente, a matriz de Routh é construída apartir do denominador da seguinte função de transferência para cadeia fechada: 𝐺 𝐾 𝐺𝐺 𝐺 (42) 𝐶 = 𝐿 + 𝑆𝑃′ 1+𝐺 𝐺 𝐺 𝐺 1+𝐺 𝐺 𝐺 𝐺 Para resolver o polinómio do denominador, considera-se 𝐺𝑣 (função transferência do atuador) e 𝐺𝑚 (função transferência do sensor) são considerados com valor 1. 𝐺𝑐 (função transferência do controlador) e 𝐺𝑝 (função de transferência da manipulada) já foram anteriormente descritas. 36  Portanto o polinómio do denominador fica como na expressão 43. 𝐾 (1 + 𝑡 𝑠) 1 𝑃(𝑠) = 1 + 𝐺 𝐺 𝐺 𝐺 = 1 + 𝑒 × 𝐾𝑐(1 + + 𝑡 𝑠) (43) (1 + 𝑡 𝑠)(1 + 𝑡 𝑠) 𝑡𝑠 A resolução deste polinómio envolve a aplicação de aproximação de Padê – Multiplicar o numerador da função de transferência da manipulada (Gm) por 1 − 𝑠 e o denominador por 1 + 𝑠. Para o sistema de tanque de nível, a resolução do denominador dá origem a um polinómio de quarto grau, em relação a s, do tipo: 𝑎 𝑠 +𝑎 𝑠 +𝑎 𝑠 +𝑎 𝑠+𝑎 =0 (43) Em que, 1 1 𝑎 =− 𝑡 𝐾𝑡 𝑡𝐾 + 𝑡 𝑡 𝑡𝜃 (44) 2 2 1 1 1 1 𝑎 =𝑡 𝐾 𝑡 𝑡 𝐾 − 𝑡 𝐾 𝑡 𝐾 𝜃− 𝐾 𝑡 𝑡 𝐾 𝜃+ 𝑡 𝑡 𝜃+𝑡 𝑡 𝑡 + 𝑡 𝑡 𝜃 (45) 2 2 2 2 1 1 1 𝑎 =𝑡 𝐾 𝑡 𝐾 − 𝑡 𝐾 𝐾 𝜃+𝐾 𝑡 𝑡 𝐾 − 𝐾 𝑡 𝐾 𝜃+𝑡 𝑡 +𝑡 𝑡 + 𝑡 𝑡 𝜃 (46) 2 2 2 1 𝑎 = 𝑡 𝐾 𝐾 +𝐾 𝑡 𝐾 − 𝐾 𝐾 𝜃+𝑡 (47) 2 𝑎 =𝐾𝐾 (48) 𝑎 𝑎 −𝑎 𝑎 𝑏 = (49) 𝑎 𝑎 𝑎 −𝑎 0 𝑏 = =𝑎 (50) 𝑎 𝑏 𝑎 −𝑏 𝑎 𝑐 = (51) 𝑏 𝑐 𝑏 −𝑏 𝑎 𝑑 = =𝑏 (52) 𝑐 Para um polinómio de quarto grau, a matriz de Routh organiza-se da seguinte forma: a4 a2 a0 a3 a1 0 b1 b2 0 c1 0 0 d1 0 0 37  Para avaliar a estabilidade do controlador obtido pelo método de Ziegler-Nicholas, esta matriz é a indicada, contudo como para a execução do modelo interno foi feita uma aproximação da função de transferência de segunda ordem com atraso para uma função de transferência de primeira ordem com atraso, a matriz a utilizar neste método é a mesma usada para o tanque de pressão. Ziegler-Nichols 1802,577 388,16517 5,593761 1147,64 62,739099 0 289,622 5,5937608 0 40,57358 0 5,593761 0 Modelo Interno 1306,552 145,2283 0 883,9777 4,138956 0 139,1107 0 0 4,138956 0 0 Assim sendo, como todos os valores da coluna mais à esquerda são positivos, ambos os controladores, cujo parâmetros foram obtidos através dos métodos para os afinar, são estáveis. 5. Aplicação de técnicas de reposta de frequências para análise de ambos os sistemas As ferramentas utilizadas para a análise da resposta em frequência dos sistemas em estudo foram o Critério de Nyquist e os Diagramas de Bode. Ambas as técnicas permitem determinar a estabilidade dos sistemas a partir do comportamento da função de transferência em open loop, oferecendo diferentes formas de visualização e interpretação. O critério de Nyquist, relaciona a função de transferência com a raiz do denominador da função, a estabilidade do sistema vai depender da localização dessas raízes no plano complexo. 38  O diagrama de Bode, por sua, vez, representa a resposta em frequências através de gráficos de magnitude e de fase. Permitindo assim uma análise empírica do ganho e da margem de fase ao longo de uma ampla faixa de frequências. 5.2 Tanque de Pressão Variável manipulada De maneira a estudar a resposta do sistema às frequências em open loop, foi feito o estudo em cadeia aberta para a variável manipulada através da análise do diagrama de Bode. Assim, este foi representado como o rácio entre a amplitude de saída da frequência aplicada e a amplitude de entrada no estado estacionário e o ângulo de fase, φ, em função das frequências aplicadas, ω. Este gráfico representa a função de transferência numa escala logarítmica. 𝑲 𝑮(𝒔) = 𝒆 𝜽𝒔 (53) 𝟏 + 𝝉𝒔 Simplificando e substituindo s por jω, temos 𝑲 (54) 𝑮𝟏 (𝒋𝝎) = 𝟏 + 𝝉𝒋𝝎 𝑮𝟐 (𝒋𝝎) = 𝒆 𝜽𝝎 (55) Sabendo que, 𝑨𝑹 = 𝑹𝟐 + 𝑰𝟐 (56) 𝟏 𝟏 (57) 𝝋 = 𝐭𝐚𝐧 𝑹 𝑲 (58) 𝑨𝑹𝟏 = 𝒆 𝝋𝟏 = −𝐭𝐚𝐧 𝟏 (𝝎𝝉) √𝝎𝟐 𝝉𝟐 + 𝟏 𝑨𝑹𝟏 = 𝟏 𝒆 𝝋𝟏 = −(𝝎𝝉) (59) Com estes parâmetros conseguimos representar graficamente os diagramas de Bode. 39  Figura 32 - Diagrama de Bode para a razão de amplitudes da variável manipulada. Figura 33 – Diagrama de Bode para o ângulo de fases da variável manipulada. Através da leitura dos gráficos conseguimos deduzir que, para frequências muito baixas, o sistema tem um ganho elevado e à medida que a frequência aumenta, o ganho decai rapidamente. Isto sugere que o sistema transmite com maior intensidade os sinais de mais baixa frequência, enquanto que sinais de alta frequência são atenuados. Já pelo diagrama de fase podemos perceber que o sistema não apresenta atraso significativo na resposta a frequências baixas. Existe um atraso significativo para sinais de alta frequência, com boa resposta a variações lentas e atenuação de variações rápidas. De seguida desenvolveu-se o diagrama de Nyquist, representado num referencial imaginário, onde AR é a distância ao centro do eixo e o ângulo de fase, φ, a distância ao eixo real positivo. 40 As partes real e imaginária foram determinadas da seguinte forma: 𝑹𝒆[𝑮(𝒋𝝎)] = 𝑨𝑹 𝐜𝐨𝐬 𝝋 (60) 𝑰𝒎[𝑮(𝒋𝝎)] = 𝑨𝑹 𝐬𝐞𝐧 𝝋 (61) Figura 34 - Diagrama de Nyquist para a variável manipulada sem atraso. Figura 35 - Diagrama de Nyquist para o tempo de atraso. 41  Figura 36 - Diagrama de Nyquist para a variável manipulada com atraso. Pela figura 34 podemos admitir que o sistema sem atraso é estável, uma vez que ao não cruzar o eixo real negativo (-1,0), não há risco imediato de instabilidade. Já a figura 35 sugere que, por ter sempre um rácio de amplitude igual a 1, o atraso não afeta diretamente o ganho, mas sim a fase, o que pode levar a instabilidade quando combinado com a dinâmica do sistema principal. Por fim, na figura 36, é possível verificar que o rácio de amplitudes estará continuamente a diminuir, a tender para 0 com o aumento das frequências e, ao mesmo tempo, a rodar em volta do centro dado que, como foi possível ver anteriormente, o ângulo de fase estará a tender para o infinito. Carga O procedimento foi o mesmo que utilizado para a variável manipulada, os resultados são apresentados a seguir. Figura 37 - Diagrama de Bode para a razão de amplitudes da carga. 42 Figura 38 - Diagrama de Bode para o ângulo de fases da carga. 43  Figura 39 – Diagrama de Nyquist para a carga sem atraso Figura 40 – Diagrama de Nyquist para o tempo de atraso. Figura 41 - Diagrama de Nyquist para a carga com atraso. As conclusões podem ser consideradas as mesmas que as anteriores. Como se trata de um sistema de primeira ordem, denota-se um decréscimo do AR e do ângulo de fase a tender para o infinito. 44 Análise da resposta do sistema às frequências em control loop (𝐺𝑂𝐿) Ziegler-Nichols – Bode e Nyquist Para realizar o estudo da resposta dos controladores às diferentes frequências aplicadas, foi necessário definir novamente AR e os ângulos de fase. Para tal, partiu-se da função de transferência de cadeira aberta: 𝐺 =𝐺 𝐺 𝐺 𝐺 (62) Tal como referido anteriormente, 𝐺𝑣 e 𝐺𝑚 apresentarão um valor de 1, portanto falta discriminar a função de transferência do controlador: 1 𝐺 (𝑠) = 𝐾 (1 + + 𝑡 𝑠) (63) 𝑡𝑠 Substituindo s por jw, 1 𝐺 (𝑗𝑤) = 𝐾 (1 + + 𝑡 𝑗𝑤) (64) 𝑡 𝑗𝑤 (𝑡 𝑤) + (1 − 𝑡 𝑡 𝑤 ) 𝐴𝑅 = 𝐾 (65) (𝑡 𝑤) (1 − 𝑡 𝑡 𝑤 ) 𝜙 = 𝑡𝑎𝑛 − (66) 𝑡𝑤 Figura 42 – Diagrama de Bode da razão de amplitude para o controlador Ziegler-Nichols, GOL e open loop. 45  Figura 43 - Diagrama de Bode do ângulo de fases para o controlador Ziegler-Nichols, GOL e open loop. Figura 44 – Diagrama de Nyquist do GOL com controlador obtido pelo Ziegler-Nichols. Tabela 15 - Frequências de corte para o controlador PID pelo Ziegler-Nichols. ωb, τI (s-1) ωb, τD (s-1) 0,57 2,29 O controlador PID aumenta o ganho para altas frequências, o que melhora a resposta rápida do sistema, mas também pode introduzir ruído. O controlador, reduz também o atraso de fase nas frequências médias, o que melhora a estabilidade e a resposta do sistema. O GOL apresenta um comportamento mais equilibrado, adequado para garantir estabilidade. 46  O controlador obtido pelo método Ziegler-Nichols estabiliza o sistema, mas com margens de estabilidade reduzidas. Isso indica que o ajuste do controlador poderia ser refinado para melhorar a robustez. Modelo Interno– Bode e Nyquist Figura 45 - Diagrama de Bode da razão de amplitude para o controlador Modelo Interno, GOL e open loop. Figura 46 - Diagrama de Bode do ângulo de fases para o controlador Modelo Interno, GOL e open loop. 47  Figura 47 - Diagrama de Nyquist do GOL com controlador obtido pelo Modelo Interno. Tabela 16 - Frequências de corte para o controlador PID pelo Modelo Interno. ωb, τI (s-1) ωb, τD (s-1) 0,21 2,28 O controlador obtido pelo Modelo Interno apresenta uma melhoria significativa em termos de robustez e estabilidade, com margens mais amplas e uma menor sensibilidade a ruídos em altas frequências. Comparado ao controlador projetado pelo método Ziegler-Nichols, o Modelo Interno oferece um desempenho mais confiável, especialmente para sistemas sujeitos a atrasos e perturbações. Estabilidade dos sistemas em Closed loop Tabela 17 – Parâmetros de estabilidade em Closed loop. Ziegler-Nichols Modelo Interno ωc (min−𝟏) 2,22 2,43 ARc 0,29 0,12 Dessa forma, segundo o critério de Bode os dois métodos apresentaram sistemas estáveis, uma vez que 𝐴𝑅𝐶 é menor que 1. Estabilidade Relativa O último método utilizado para estudar a estabilidade do sistema foi a estabilidade relativa, a qual é realizada através da margem de ganho (GM) e a margem de fase (PM) 48  1 (67) 𝐺𝑀 = 𝐴𝑅 𝑃𝑀 = 180 + 𝜙 (68) Onde 𝐴𝑅𝐶 é a amplitude crítica e 𝜙 é ângulo de fase obtido para a frequência critica do ganho ωg, para o qual o 𝐴𝑅GO𝐿 é 1. Tabela 18 - Resultados do estudo da estabilidade relativa para ambos os controladores. Ziegler-Nichols Modelo Interno ωc (s-1) 2,44 2,56 Arc 0,29 0,12 Ωg (s-1) 0,44 0,39 𝝓g -135,1º -102,5º GM 3,44 8,33 PM 44,9º 77,5º Embora o Modelo Interno ofereça maior estabilidade e robustez, a margem de ganho apresentada supera o intervalo de valores recomendados (1,7 - 4), assim indica que este pode ser mais lento na resposta do sistema, apesar de apresentar outros parâmetros de estabilidade superior ao do Ziegler-Nichols que apesar de oferecer uma resposta mais rápida, apresenta menor robustez. 5.3 Tanque de Nível Variável Manipulada Tal como foi realizado para o tanque de pressão, primeiramente foi realizado o estudo em cadeia aberta para a variável manipulada. Para tal, foi designado o diagrama de bode tendo por base a função de transferência de segunda ordem que representa o sistema. Para a construção do diagrama, é necessário dividir a função de transferência nos seus parâmetros, substituindo s por jω, para que se possa estabelecer a parte real e a parte imaginária de cada parâmetro. 𝐺1(𝑗𝜔) = 𝐾 (68) 𝐺2 (𝑗𝜔) = 𝑡 𝑗𝜔 (69) 1 𝐺3(𝑗𝜔) = (70) 1 + 𝑡 𝑗𝑤 49  1 𝐺4(𝑗𝜔) = (71) 1 + 𝑡 𝑗𝑤 𝐺5(𝑗𝜔) = 𝑒 (72) Tendo por base que, 𝐴𝑅 = 𝑅 +𝐼 (73) 𝐼 𝜙 = 𝑡𝑎𝑛 (74) 𝑅 Podemos calcular AR - o rácio entre a amplitude de saída consequente da frequência aplicada e a amplitude de entrada no estado estacionário – e o ângulo de fase (φ) para cada um dos parâmetros. 𝐴𝑅 = 𝐾 𝜙 =0 (75) 𝐴𝑅 = 𝑤 𝑡 +1 𝜙 = 𝑡𝑎𝑛 (𝑤𝑡 ) (76) 𝐴𝑅 = 𝜙 = − 𝑡𝑎𝑛 (𝑤𝑡 ) (77) 𝐴𝑅 = 𝜙 = − 𝑡𝑎𝑛 (𝑤𝑡 ) (78) 𝐴𝑅 = 1 𝜙 = −(𝑤𝜃) (79) 100000000 1000000 10000 100 AR 1 0,01 0,0001 0,000001 1E-08 1E-10 0,001 0,01 0,1 1 10 100 1000 10000 100000 AR1 w AR2 AR3 AR4 AR5 AR(GOL) com atraso AR (GOL) sem atraso Figura 48 - Diagrama de Bode relativo às razões de amplitude para a variável manipulada 50  200 100 0 φ -100 -200 -300 -400 0,001 0,01 0,1 1 10 100 1000 10000 100000 φ1 w φ2 φ3 φ4 φ5 φ (GOL) com atraso φ (GOL) sem atraso Figura 49 - Diagrama de Bode para os diversos ângulos de fase da variável manipulada Pela análise do diagrama de Bode, o aumento da frequência provoca alterações na maior parte dos declives das curvas de razões de amplitude: as curvas de t1 e t2 (AR3 e AR4) começam a apresentar declives negativos, enquanto que a curva de ta apresenta um valor positivo. Como seria expectável tanto a curva do atraso e da constante do ganho (AR5 e AR1) mantém-se constantes. A análise do gráfico permite também concluir que as curvas do ta e t1 possam ser consideradas simétricas, significando que estes dois parâmetros se possam anular. Com esta suposição, é possível assumir que a cadeia aberta apresentará um comportamento de um sistema de primeira ordem, com declive negativo. O diagrama de ângulo de fase obtido indica que com o aumento da frequência, ta descreve um declive positivo até tender para a zona dos 90º e em concordância, com o diagrama de Bode, τ1 e τ2 apresentam declives negativos até tenderem para a zona dos -90º. Posto isto, se a curva do K não apresenta declive, a curva para a cadeia aberta sem atraso terá declive negativo. Em relação à curva do atraso, esta contém um declive negativo e que tende para menos infinito. Como tal, é percetível que a curva da cadeia aberta com atraso apresente um declive negativo e que contenha valores cada vez mais negativos. Tal como no tanque de pressão os pontos de interseção dos declives são caracterizados pelas chamadas frequências de corte, 𝜔𝑏, neste caso, como é um sistema de segunda ordem com dinâmica de numerador, existe a influência de 3 tempos característicos, pelo que existirão 3 frequências de corte. 51  Tabela 19 - Frequências de corte para a variável manipulada. 𝜔𝑏,𝜏1(𝑚𝑖𝑛−1 ) 𝜔𝑏,𝜏a(𝑚𝑖𝑛−1 ) 𝜔𝑏,𝜏2(𝑚𝑖𝑛−1 ) 0,28 0,24 0,04 Para o sistema de tanque de nível, a realização do diagrama de Nyquist designou-se para três situações: para a variável manipulada com atraso, para o atraso e para a variável manipulada sem atraso. A análise do diagrama de Nyquist para a variável manipulada com atraso segui as seguintes equações: 𝑅𝑒 𝐺(𝑗𝑤) = 𝐴𝑅 ( ) cos (𝜙 ( )) (80) 𝐼𝑚 𝐺(𝑗𝑤) = 𝐴𝑅 ( ) sin (𝜙 ( )) (81) 0,1 0 -0,1 -0,2 lm -0,3 -0,4 -0,5 -0,6 -0,2 0 0,2 0,4 0,6 0,8 1 Re Figura 50 - Diagrama de Nyquist para a variável manipulada com atraso Para o atraso, 𝑅𝑒[𝐺(𝑗𝜔)] = 𝑅𝑒(𝑒 ) = 𝑐𝑜𝑠 (𝜃𝑤) (82) 𝐼𝑚[𝐺(𝑗𝜔)] = 𝐼𝑚(𝑒 ) = 𝑠𝑒𝑛 (𝜃𝑤) (83) 1,5 1 0,5 0 lm -1,5 -1 -0,5 0 0,5 1 1,5