PID digital – Método de cálculo numérico

Caro Colega,

controle PID é uma das formas mais conhecidas e tradicionais de se implementar controle de processos. O desafio é embutir esse controle em microcontroladores e de baixo custo e com recursos de hardware limitados. No artigo técnico Controlador PID digital: Uma modelagem prática para microcontroladores – Parte 1, de Felipe Neves,  publicado originalmente no site Embarcadosé  apresentado com detalhes o controlador PID, sua arquitetura típica, o desenvolvimento das equações para o controle e alguns exemplos de como se pode desenvolver um controle PID. Nesse artigo, também publicado originalmente no site Embarcados,  serão retomados alguns conceitos sobre esse tipo de controle e serão apresentadas formas de cálculo numérico para implementar a parte diferencial e a integral do controlador, de forma a viabilizar esses cálculos em microcontroladores. No final são apresentados códigos que você pode executar diretamente no Scilab para experimentar e explorar os detalhes dos resultados obtidos. Divirta-se!

Abraço,

 

Henrique

consulte sempre um engenheiro eletrônico

======================================================

PID_Capa

Introdução

 

No artigo técnico Controlador PID digital: Uma modelagem prática para microcontroladores – Parte 1é  apresentado com detalhes o controlador PID, sua arquitetura típica, o desenvolvimento das equações para o controle e alguns exemplos de como se pode desenvolver um controle PID. Nesse artigo serão retomados alguns conceitos sobre esse tipo de controle e serão apresentadas formas de cálculo numérico para implementar a parte diferencial e a integral do controlador. Essa forma alternativa de cálculo casa muito bem com a aplicação em microcontroladores ou processadores de 8 bits ou outros relativamente lentos. Os resultados que serão apresentados a seguir, foram largamente aplicados nos primeiros controladores PID digitais e certamente ainda são úteis nos dias de hoje, por sua simplicidade e velocidade de processamento.

Arquitetura típica de um controlador PID

 

Um sistema de controle típico de PID em malha fechada tem o aspecto mostrado na Figura 1. Nessa figura pode-se observar uma entrada, referente a um sinal de referência ou setpoint, um bloco onde essa referência é subtraída de um sinal proveniente do processo, em geral medido de alguma forma por um sensor, transdutor etc. Essa subtração produz um sinal de erro que entra no controlador PID propriamente dito, e a saída do controlador, que atua sobre o processo que está sendo controlado, é aplicada para conduzir o processo ao estado desejado.

.

PID_1a

 Figura 1: Controle PID típico

 

Na prática, o controlador PID digital destacado na Figura 1 em vermelho, pode ser representado como mostrado na Figura 2. O erro é amostrado, produzindo assim um conjunto de amostras discretizadas do sinal de erro, que é conduzido a um conversor analógico / digital (A/D), depois são calculados os valores da derivada, da integral e do valor proporcional, em seguida esses valores calculados são somados e conduzidos a um conversor digital / analógico (D/A).  

PID_2a

Figura 2: Controlador PID em blocos

 

O desenho da Figura 2 representa o fluxo do sinal de erro para o controlador PID de forma a simplificar o desenvolvimento dos cálculos que serão apresentados a seguir. É evidente que num controlador real, os sinais de SP (setpoint) e PV (present value) também são digitalizados e subtraídos para que se obtenha o erro. O importante é fixar o conceito de que o erro deve ser entregue ao controlador numa frequência ou taxa fixa. Isso é requisito para o desenvolvimento que será realizado a seguir. Para que o sistema funcione bem, também é necessário que a atualização do sinal de saída pelo conversor D/A seja feito na mesma frequência do amostrador de entrada. Na Figura 3 é mostrado um sinal analógico como exemplo de sinal de erro e são destacadas as amostras realizadas num período fixo (Δ).

PID_3

Figura 3: Exemplo de sinal de erro e de sua amostragem

A seguir será desenvolvido um pouco de teoria, o suficiente para que se possa saber de onde foram tiradas as equações que resultarão desse desenvolvimento. Se você preferir, pode consultar direto os quadros resumo que apresentam as equações resultantes desses cálculos.

 

Séries de Taylor

 

Brook Taylor elaborou uma teoria a respeito de funções matemáticas em 1715, conhecida como teorema de Taylor ou séries de Taylor. Esse teorema estabelece que: sendo f  uma função  infinitamente diferenciável, ela pode ser representada pela seguinte fórmula:

Equação_1a(1)

onde:

  • fx0a  representa a k-ésima derivada da função avaliada no ponto x0;
  • X0             é um ponto arbitrário de referência;
  • k              é o fatorial de k.

Se a somatória acima convergir, a função pode ser representada por uma série finita de termos. Como por exemplo a expansão da série (1):

seriea(2)

 

Obs: Quando x0 = 0, a série (2) é conhecida como série de McLaurin.

A seguir serão utilizadas as séries de Taylor para deduzir formas simples para o cálculo da primeira derivada de uma função e a integral dessa função.

 

Derivação Numérica

 

Seja definida uma função f(x), representada por uma série de Taylor, tal como mostrada em (1), onde x = i é uma amostra arbitrária de índice i tal como mostrado na Figura 3. O valor de f(x) será batizado de fi. Seja f(i-1) o valor dessa mesma função calculado no índice (i-1). A distância entre essas amostras é de Δ, correspondente a uma unidade de tempo relacionada com a frequência ou taxa de amostragem. Daí obtemos:

Equação_fi(3)

Para avaliar a função em (i-1), vamos substituir i por (i-Δ). Realizando a substituição e rearranjando a equação obtida, chega-se na equação (4).

Equação_fi_1(4)

Escrevendo essa somatória na forma expandida, obtém-se:

Equação_fi_1a (5)

Se desprezarmos os termos de segunda ordem e superiores e rearranjando, a equação fica assim:

Equação_2pontos(6),  e o erro = erro_2pontos

 

A expressão (6) é conhecida como a derivada com 2 pontos.

 

É evidente que nesse cálculo, o erro é relativamente grande. Pode-se estender o raciocínio desenvolvido acima para um cálculo mais preciso, que pode ser obtido com 3, 4 ou mais amostras do sinal.

 

Se substituirmos i por (i- 2*Δ) e for desenvolvido o mesmo raciocínio que foi utilizado para i-1 poderemos levantar a equação de f(i-2) em função de fi.

Equação_fi_2(7)

Usando a equação (5) e rearranjando, obtém-se:

 

derivada_3pontos(8), onde o erro = erro_3pontos

 

Essa é a equação de cálculo de derivada com 3 pontos. De maneira análoga pode-se desenvolver a derivada de 4 pontos, mostrada na equação (9).

 

derivada_4pontos(9), onde o erro = erro_4pontos

 

A seguir será apresentado um quadro resumo com as fórmulas desenvolvidas até aqui.

 

Quadro resumo da diferenciação digital

 

Aqui será adotada uma nova notação, mais frequentemente utilizada em processamento digital de sinais (DSP) e em controle digital. Isso facilitará o entendimento de quem está habituado a essa notação e introduz a notação para quem ainda não a conhece. A notação é assim:

  • y = saída (nesse caso a derivada da função);
  • x = amostra atual da função;
  • z-n = n-ésima amostra anterior;
  • Δ = intervalo de tempo entre amostras.

Resumo das fórmulas referentes à diferenciação digital:

 

=====================================================================

  • Derivada calculada com 2 pontos:

der_2p_z(10)

—————————————————————————————————————————

  • Derivada calculada com 3 pontos:

der_3p_z(11)

————————————————————————————————————————–

  • Derivada calculada com 4 pontos:

der_4p_z(12)

===================================================================

 

Os resultados obtidos nas equações acima são facilmente implementáveis em microprocessadores.

 

 

Integração Numérica

 

Para desenvolver as equações referentes à integração numérica, o raciocínio é o mesmo já desenvolvido para a diferenciação. Voltando à Figura 3, onde temos representada uma forma de onda arbitrária amostrada em intervalos regulares de Δ segundos, pode-se definir que:

Integral_0(13)

Seja o segundo termo da equação (13):

integral0a(14)

Desenvolvendo a série para Ii-1 se obtém uma equação semelhante à equação (5):

integral1(15)

Rearranjando, obtém-se:

Integral2(16)

Desenvolvendo:

integral3(17)

Substituindo:

integral4(18)

desprezando os termos de segunda ordem ou superiores chega-se à equação para o cálculo da integral com 2 pontos.

 

integral5(19), onde o erro = erro_integral

 

A fórmula (19) é conhecida como “regra do trapézio”.

 

O cálculo da integral pode ser estendido para mais pontos, para aumentar a precisão do resultado. A seguir, a equação para o cálculo com 3 pontos (20) e 4 pontos (21).

 

integral6(20), onde o erro =erro3i

.

integral7(21), onde o erro = erro4i

 

Quadro resumo da integração digital

Aqui também será adotada a nova notação, mais frequentemente utilizada em processamento digital de sinais (DSP) e em controle digital.  A notação é assim:

  • y = saída (nesse caso a integral da função);
  • x = amostra atual da função;
  • z-n = n-ésima amostra anterior;
  • Δ = intervalo de tempo entre amostras.

Resumo das fórmulas referentes à integração digital:

 

=====================================================================

  • Integral calculada com 2 pontos:

Integralz2p(22)

 

—————————————————————————————————————————

  • Integral calculada com 3 pontos:

integral3pz(23)

 

————————————————————————————————————————–

  • Integral calculada com 4 pontos:

integral4pz(24)

===================================================================

 

Observe que as fórmulas acima calculam a integral incremental entre uma amostra e a anterior. Para que seja realizado um cálculo de integral de fato é necessário acumular os resultados.

 

Obtidos os resultados  acima, é interessante verificar se eles realmente funcionam. Para isso eles serão simulados num programa para cálculos matemáticos conhecido por Scilab. Serão mostrados os respectivos códigos e os resultados. Para que você possa testar e explorar na sua própria máquina será necessário instalar o Scilab nela. Baixe o programa aqui: SCILAB vers. 5.5.0. O programa é gratuito e tem versões para Windows, GNU/Linux e Mac OS X.

Inicie o programa e digite o código a seguir. Se você preferir, você pode marcar o código, copiá-lo (Ctrl-C) e colá-lo no Scilab (Ctrl-V). Para iniciar as operações basta teclar o “Enter”.

 

Cálculo de Integrais numéricas utilizando o Scilab

 

// Inicializar as variáveis e a forma de onda para os testes

nCiclos = 10; // Define o núemro de ciclos da forma de onda
nAmostras = 40; // Define o número de amostras por ciclo

// Gera o vetor de referência para a produção da forma de onda

t = (1:2*%pi/nAmostras:nCiclos * 2 * %pi + (1-2*%pi/nAmostras)); 
nvErro = squarewave(t); // Gera uma onda quadrada com 400 amostras

// Gera o gráfico para essa forma de onda

a0 = figure(0);
a0.background = -2;

plot(nvErro);

// Gera o grid e as legendas

xtitle('Forma de onda do Erro para testes', 'Indice das amostras', 'Amplitude');
set(gca(), "data_bounds", [0,400,-1.2,1.2]);
set(gca(), "grid", [1,1]);

 

 

Os comandos acima geram um gráfico parecido com o apresentado na Figura 4.

PID_Fig_1a

Figura 4: Gráfico do forma de onda quadrada gerada no Scilab

 

Como não há formas de onda triangulares prontas no Scilab, será utilizada a forma de onda quadrada, calculada sua integral para formar a onda triangular e calculada a diferencial da onda triangular para formar a onda quadrada novamente. Digite o seguinte código no Scilab (ou copie e cole):

 

// Inicializa as variáveis para o cálculo da integral com 2 pontos

Int2p = zeros(1, size(nvErro,2));
Int2pAcum = 0;
Z1 = 0;
Delta = 1;

// Calcula a Integral com 2 pontos


for nIndice = 1:size(nvErro,2)
    Int2p(nIndice) = ((nvErro(nIndice) + Z1)* Delta/2) + Int2pAcum;
    Int2pAcum = Int2p(nIndice);
    Z1 = nvErro(nIndice);
end;

// Apresenta a integral no gráfico

a1 = figure(1);
a1.background = -2;

plot(Int2p, 'm');

// Gera o grid e as legendas

set(gca(), "background", [-2]); 
xtitle('Integral da forma de onda calculada com 2 pontos', 'Indice das amostras', 'Amplitude'); 
set(gca(), "data_bounds", [0,400,-10,20]); 
set(gca(), "grid", [1,1]);

//==========================================================================

// Inicializa as variáveis para o cálculo da integral com 3 pontos

Int3p = zeros(1, size(nvErro,2));
Int3pAcum = 0;
Z1 = 0;
Z2 = 0;
Delta = 1;

// Calcula a Integral com 3 pontos

for nIndice = 1:size(nvErro,2)
     Int3p(nIndice) = ((5*nvErro(nIndice) + 8*Z1- Z2)* Delta/12) + Int3pAcum;
     Int3pAcum = Int3p(nIndice);
     Z2 = Z1;
     Z1 = nvErro(nIndice);
end;

// Gera gráfico

a2=figure(2);
a2.background = -2;

plot(Int3p, 'b');

// Gera o grid e as legendas

set(gca(), "background", [-2]); 
xtitle('Integral da forma de onda calculada com 3 pontos', 'Indice das amostras', 'Amplitude'); 
set(gca(), "data_bounds", [0,400,-10,20]); 
set(gca(), "grid", [1,1]);

//========================================================================================================= // Inicializa as variáveis para o cálculo da integral com 4 ponto

Int4p = zeros(1, size(nvErro,2));
Int4pAcum = 0; 
Z1 = 0; 
Z2 = 0; 
Z3 = 0; 
Delta = 1; 

// Calcula a Integral com 4 pontos 

for nIndice = 1:size(nvErro,2) 
    Int4p(nIndice) = ((9*nvErro(nIndice) + 19*Z1- 5*Z2 + Z3)* Delta/24) + Int4pAcum;
    Int4pAcum = Int4p(nIndice); 
    Z3 = Z2; 
    Z2 = Z1; 
    Z1 = nvErro(nIndice);
end; 

// Gera gráfico figure(3);

a3 = figure(3);
a3.background = -2;

plot(Int4p, 'r');

// Gera o grid e as legendas

set(gca(), "background", [-2]); 
xtitle('Integral da forma de onda calculada com 4 pontos', 'Indice das amostras', 'Amplitude'); 
set(gca(), "data_bounds", [0,400,-10,20]); 
set(gca(), "grid", [1,1]); 

 

Os gráficos obtidos são como o mostrado na Figura 5.

Integral_1a

Figura 5: Gráfico da integral calculada com 3 pontos

 

Os três gráficos obtidos parecem ser iguais, por motivo da forma de onda quadrada ser muito simples. Se olharmos com detalhes os “bicos” das formas de onda ampliadas, poderemos observar as diferenças (Figura 6). Quanto maior o número de pontos usados para o cálculo, menor a diferença se comparada com a onda triangular ideal.

Integral_detalhe

Figura 6: Gráfico das integrais calculadas sobrepostas

 

Cálculo de Diferencias numéricas utilizando o Scilab

 

Em seguida serão mostradas a utilização das diferenciais calculadas com 2, 3 ou 4 pontos. A onda utilizada para essas operações é a onda triangular resultante da  integral calculada com 4 pontos.

 

// Calculo das diferencias numéricas utilizando-se 2 pontos

// Fecha os gráficos

close(a0);  
close(a1);
close(a2);
close(a3);


Diff2p = zeros(1, size(Int4p,2));
Z1 = 0;
Delta = 1;

// Calcula a Diferencial com 2 pontos

for nIndice = 1:size(nvErro,2)
    Diff2p(nIndice) = ((Int4p(nIndice) - Z1)* 1/Delta);
    Z1 = Int4p(nIndice);
end;

// Apresenta a diferencial no gráfico

a0 = figure(0);
a0.background = -2;

plot(Diff2p, 'm');

// Gera o grid e as legendas

set(gca(), "background", [-2]); 
xtitle('Diferencial da forma de onda calculada com 2 pontos', 'Indice das amostras', 'Amplitude'); 
set(gca(), "data_bounds", [0,400,-1.5,1.5]); 
set(gca(), "grid", [1,1]);

//==========================================================================

// Calculo das diferencias numéricas utilizando-se 3 pontos

Diff3p = zeros(1, size(Int4p,2));
Z1 = 0;
Z2 = 0;
Delta = 1;

// Calcula a Diferencial com 3 pontos

for nIndice = 1:size(nvErro,2)
    Diff3p(nIndice) = ((3*Int4p(nIndice) - 4*Z1 + Z2)* 1/(2*Delta));
    Z2 = Z1;
    Z1 = Int4p(nIndice);
end;

// Apresenta a deferencial no gráfico

a1 = figure(1);
a1.background = -2;

plot(Diff3p, 'b');

// Gera o grid e as legendas

set(gca(), "background", [-2]); 
xtitle('Diferencial da forma de onda calculada com 3 pontos', 'Indice das amostras', 'Amplitude'); 
set(gca(), "data_bounds", [0,400,-2.5,2.5]); 
set(gca(), "grid", [1,1]);

//==========================================================================

// Calculo das diferencias numéricas utilizando-se 4 pontos

Diff4p = zeros(1, size(Int4p,2));
Z1 = 0;
Z2 = 0;
Z3 = 0;
Delta = 1;

// Calcula a Diferencial com 4 pontos

for nIndice = 1:size(nvErro,2)
    Diff4p(nIndice) = ((11*Int4p(nIndice) - 18*Z1 + 9*Z2 - 2*Z3)* 1/(6*Delta));
    Z3 = Z2;
    Z2 = Z1;
    Z1 = Int4p(nIndice);
end;

// Apresenta a deferencial no gráfico

a2 = figure(2);
a2.background = -2;

plot(Diff4p, 'r');

// Gera o grid e as legendas

set(gca(), "background", [-2]); 
xtitle('Diferencial da forma de onda calculada com 4 pontos', 'Indice das amostras', 'Amplitude'); 
set(gca(), "data_bounds", [0,400,-2.5,2.5]); 
set(gca(), "grid", [1,1]);

 

Os gráficos obtidos são semelhantes ao da Figura 7.

 Diferencial_1b

Figura 7: Diferencial da onda triangular calculada com 4 pontos

 

Na Figura 8 pode-se observar as três diferenciais sobrepostas. Observe que a diferencial com mais pontos tem um transitório maior do que as outras, o que está coerente. Aqui encerram-se as simulações para demonstrar o funcionamento dos métodos numéricos.

É recomendável que o leitor faça experimentos com o Scilab, crie outras formas de onda, utilize formas de onda reais, pois o Scilab permite que se carregue dados de arquivos, e se familiarize com os recursos do Scilab, se acaso ainda não os conhece. O Scilab é uma ferramenta poderosa, tem recursos para cálculos de problemas de controle, processamento digital de sinais entre muitos outros.

Diferencial_sobrep

Figura 8: Detalhe da sobreposição das diferenciais calculadas

Nesse artigo técnico foram desenvolvidos métodos numéricos para o cálculo de integrais e diferenciais, típicas para uso em controle digital por PID. Esses métodos foram simulados no Scilab e foi mostrado o funcionamento correto das diversas aproximações e uma melhor precisão à medida que se utiliza mais pontos para o cálculo dos resultados. As equações desenvolvidas são muito adequadas para a utilização em microprocessadores pela sua simplicidade. Experimente usá-las.

Licença Creative Commons

Esta obra, “ PID digital – Método de cálculo numérico“, de Henrique Frank W. Puhlmann, foi licenciada sob uma Licença Creative Commons Atribuição-Não Comercial-CompartilhaIgual 3.0 Não Adaptada.

Um comentário sobre “PID digital – Método de cálculo numérico

Deixe um comentário

Preencha os seus dados abaixo ou clique em um ícone para log in:

Logotipo do WordPress.com

Você está comentando utilizando sua conta WordPress.com. Sair /  Alterar )

Foto do Google

Você está comentando utilizando sua conta Google. Sair /  Alterar )

Imagem do Twitter

Você está comentando utilizando sua conta Twitter. Sair /  Alterar )

Foto do Facebook

Você está comentando utilizando sua conta Facebook. Sair /  Alterar )

Conectando a %s