P ROGRAMA DE P ÓS G RADUAÇÃO EM E NGENHARIA
ANÁLISE NO DOMÍNIO DA FREQUÊNCIA PARA
AUXÍLIO NA SELEÇÃO DE ESTRUTURA DE MODELOS
NARX POLINOMIAIS
Felipe de Brito Freitas
Dissertação submetida à banca examinadora designada pelo
Colegiado do Programa de Pós-Graduação em Engenharia do
Centro Universitário do Leste de Minas Gerais, como parte
dos requisitos necessários à obtenção do grau de Mestre em
Engenharia Industrial.
Área de Concentração: Processos Industriais
Orientador:
Prof. Marcelo Vieira Corrêa, Dr. - PPGE/Unileste-MG
Coronel Fabriciano, março de 2011
P ROGRAMA DE P ÓS G RADUAÇÃO EM E NGENHARIA
ANÁLISE NO DOMÍNIO DA FREQUÊNCIA PARA
AUXÍLIO NA SELEÇÃO DE ESTRUTURA DE MODELOS
NARX POLINOMIAIS
Felipe de Brito Freitas
Banca:
Prof. Marcelo Vieira Corrêa, Dr. - PPGE/Unileste-MG - Orientador.
Prof. Mara Cristina da Silveira Coellho, Dr. - UNIFEI.
Prof. Roselito de Albuquerque Teixeira, Dr. - PPGE/Unileste-MG.
À Deus,
à meus pais Manoel Ricardo Pacheco e Eula e
à todos que torceram por mim.
Agradecimentos
Agradeço a Deus por permitir que eu realize este trabalho dando-me condições físicas, emocionais e espirituais durante a minha vida para a realização deste trabalho.
Aos meus pais, Manoel Ricardo e Eula Elizabeth pelo apoio incondicional que me concedem a
cada momento provendo a mim todos os recursos para que eu supere cada etapa em busca da
realização dos meus sonhos.
Agradeço aos colegas do MOCP, em especial aqueles que conviveram comigo durante os dois
anos de trabalho. Citando em especial os bolsistas de iniciação científ ca Igor, Renan, Marcelo
e Natália, Tallys e também os colegas de Mestrado, Raquel, Renata, Marcus, Viviane, Gabriela,
Marisa, Felipe e Maressa que apoiaram e criticaram nos momentos certos.
À todos os professores do Mestrado que tiveram a paciência do ensinar e colaboraram de forma
imensurável para a realização deste trabalho. Em especial aqueles que tive o prazer de conviver
durante o curso de suas disciplinas: Mário Godinho, Roselito, Adelaide, Marcelo e Dair. Àqueles professores do programa que não tive a oportunidade de conviver em sala de aula mas que
sempre que precisei deles estavam de prontidão.
Às secretárias do programa que sempre pacientemente nos atendiam, esclareciam nossas dúvidas e alertavam sobre as informações importantes durante o período de curso.
Aos funcionários do UnilesteMG que sempre estiveram prestando os seus serviços de forma a
tornar o nosso trabalho menos árduo.
Ao professor Roselito de Albuquerque Teixeira pelas orientações e pela competência ao resolver
os problemas do programa na posição de coordenador do curso.
Aos familiares, amigos, vizinhos e todos aqueles que estiveram ao meu redor durante este
período sempre com palavras de apoio e incentivo.
Agradeço ao professor Marcelo Vieira Corrêa, meu orientador, pela dedicação e competência,
pelos ensinos e pela compreensão.
vi
E f nalmente, um agradecimento especial à CAPES e ao UnilesteMG pelo f nanciamento e por
depositar em mim a conf ança na minha capacidade para a realização deste trabalho.
viii
"O ser humano vivência a si mesmo, seus pensamentos como algo separado do resto
do universo - numa espécie de ilusão de ótica de sua consciência. E essa ilusão é
uma espécie de prisão que nos restringe a nossos desejos pessoais, conceitos e ao
afeto por pessoas mais próximas. Nossa principal tarefa é a de nos livrarmos dessa
prisão, ampliando o nosso círculo de compaixão, para que ele abranja todos os seres
vivos e toda a natureza em sua beleza. Ninguém conseguirá alcançar completamente
esse objetivo, mas lutar pela sua realização já é por si só parte de nossa liberação e
o alicerce de nossa segurança interior."
Albert Einstein
Resumo
Diante das dif culdades apresentadas na etapa de seleção de estrutura em identif cação de sistemas, este trabalho investiga qual a relação entre a resposta em frequência de sistemas não
lineares e a escolha de estruturas e de grau de não linearidade de modelos que o representem.
Contribuiu-se com a proposição de uma metodologia que estime estruturas adequadas à uma
representação polinomial a partir da resposta em frequência dos sinais de entrada e saída. O
procedimento baseia-se na resposta em frequência da saída gerada pelo sistema a ser modelado
excitado por um sinal de entrada de frequência única. Diante desse resultado, sugere-se o grau
de não linearidade para um modelo e determina-se os máximos atrasos de entrada e saída. O
critério ERR é utilizado para ordenar os agrupamentos de termos em ordem de relevância na
explicação do sistema. De um em um, agrupamentos são sugeridos ao modelo e verif ca-se sua
resposta em frequência respeitando a hierarquia imposta pelo ERR. Assim que há uma aproximação satisfatória da resposta em frequência do modelo em relação à resposta em frequência
do sistema é feito o corte na estrutura. O modelo é então validado. O procedimento pode tanto
ser utilizado como ferramenta única na etapa de escolha de estrutura como ferramenta auxiliar
à procedimentos clássicos existentes na literatura como o critério de ordenação ERR e o critério
de corte de AKAIKE. O procedimento proposto foi aplicado a dois sistemas simulados e a dois
sistemas reais. O primeiro sistema simulado apresenta uma não linearidade branda. O segundo
trata-se de uma planta de pH simulada com uma não linearidade severa. Observou-se que o
procedimento proposto apresentou resultados satisfatórios tanto em relação ao desempenho na
simulação dos modelos obtidos quanto à obtenção de estruturas mais compactas sem grandes
perdas de desempenho. Os sistemas reais utilizados tratam de um forno de bancada e o segundo
trata-se de um protótipo de helicóptero de dois graus de liberdade Quanserr . Dados foram coletados e o método foi aplicado obtendo resultados satisfatórios. Foram obtidas estruturas mais
compactas com pouca redução no desempenho dos modelos. Comparações foram realizadas
entre os modelos obtidos pelas técnicas clássicas de seleção de estrutura e o procedimento proposto verif cando-se ef cácia do método.
Abstract
Given the diff culties presented in the structure selection stage in system identif cation, this
work investigates what is the relationship between the nonlinear systems frequency response
and the choice of structure and the degree of nonlinearity of models that represent them. The
development of a methodology to estimate the structures is needed for a polynomial representation from the frequency response of input and output. The procedure is based on the frequency
response of the output generated by the system being modeled excited by an input signal of
single frequency. Given this result, it is suggested that the degree of nonlinearity for a model
and determines the maximum delay of input and output. The ERR criterion is used to sort
groups of terms in order of relevance in explaining the system. By one to one, groups of terms
are suggested to the model and verif es its frequency response respecting the hierarchy imposed
by the ERR. So there is a satisfactory approximation of the frequency response of the model
relative to the frequency response of the system is cut in the structure. The model is then validated. The procedure can be used both as a unique tool in the stage of choice of structure as
an auxiliary tool to conventional procedures in the literature as the ERR and Akaike criteria.
The proposed procedure was applied to two simulated systems and two real systems. The f rst
simulated system has a soft nonlinearity. The second is a plant with a pH simulated severe nonlinearity. It was observed that the proposed procedure with satisfactory results both in relation
to the performance obtained in simulation models as to obtain more compact structures without
signif cant loss of performance. The real systems used to treat a countertop oven and the second
it is a prototype helicopter to two degrees of freedom Quanser r . Data were collected and the
method was applied satisfactory outcomes. More compact structures were obtained with little
reduction in performance of the models. Comparisons were made between models obtained by
classical techniques of selection procedure proposed structure and verifying the effectiveness
of the method.
Sumário
Resumo
xi
Abstract
xiii
Lista de Figuras
xvii
Lista de Tabelas
xxi
Lista de Símbolos
1 Introdução
xxiv
1
1.1
Motivação . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
1.2
Objetivos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
1.3
Organização do texto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
2 Fundamentação Teórica
2.1
Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
9
9
xiv
2.2
Identif cação de sistemas . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.2.1
2.3
2.4
9
Etapas da Identif cação . . . . . . . . . . . . . . . . . . . . . . . . . .
10
Transformada de Fourier . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
18
2.3.1
As Transformadas de Fourier no tempo discreto . . . . . . . . . . . . .
19
2.3.2
Efeitos da amostragem de sinais . . . . . . . . . . . . . . . . . . . . .
21
Considerações Finais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
23
3 Seleção de estrutura através da análise de resposta em frequência
25
3.1
Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
25
3.2
Análise matemática da resposta em frequência de modelos NARX polinomiais .
26
3.2.1
Agrupamento de termos puros de entrada . . . . . . . . . . . . . . . .
26
3.2.2
Agrupamento de termos cruzados . . . . . . . . . . . . . . . . . . . .
32
3.2.3
Agrupamento de termos puros de saída . . . . . . . . . . . . . . . . .
35
3.2.4
Considerações sobre a análise matemática no domínio da frequência de
3.3
modelos NARX polinomiais . . . . . . . . . . . . . . . . . . . . . . .
39
Metodologia aplicada . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
41
3.3.1
Sugestão de grau de não linearidade ` a partir da análise no domínio da
frequência . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.3.2
3.4
41
Seleção de estrutura a partir da resposta em frequência da saída do sistema 43
Comentários Finais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
47
xv
4 Estudo de casos simulados
49
4.1
Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
49
4.2
Sistema simulado 1 - CTV . . . . . . . . . . . . . . . . . . . . . . . . . . . .
49
4.2.1
Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 1 utilizando o método clássico e a sugestão de grau de não
linearidade proposta . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2.2
4.3
53
Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 1 utilizando a resposta em frequência . . . . . . . . . . . . . .
57
Sistema simulado 2 - A planta de neutralização de pH simulada . . . . . . . . .
65
4.3.1
Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 2 utilizando o método clássico e a sugestão de grau de não
linearidade proposta . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.3.2
Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 2 utilizando a resposta em frequência . . . . . . . . . . . . . .
73
Validação estática dos modelos . . . . . . . . . . . . . . . . . . . . . .
75
Comentários Finais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
77
4.3.3
4.4
68
5 Estudo de casos reais
79
5.1
Introdução . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
79
5.2
Sistema Térmico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
79
5.2.1
80
Aquisição de dados . . . . . . . . . . . . . . . . . . . . . . . . . . . .
xvi
5.2.2
Obtenção de estrutura de modelo polinomial NARX para o sistema térmico utilizando o método clássico e a sugestão de grau de não linearidade proposta . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.2.3
Obtenção de estrutura de modelo NARMAX polinomial para o sistema
térmico utilizando a resposta em frequência . . . . . . . . . . . . . . .
86
Validação estática dos modelos . . . . . . . . . . . . . . . . . . . . . .
87
O helicóptero de 2 graus de liberdade . . . . . . . . . . . . . . . . . . . . . .
88
5.3.1
Aquisição de dados . . . . . . . . . . . . . . . . . . . . . . . . . . . .
90
5.3.2
Obtenção de estrutura de modelo NARMAX polinomial para o he-
5.2.4
5.3
82
licóptero utilizando o método clássico e a sugestão de grau de não linearidade proposta . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
5.3.3
5.4
Obtenção de estrutura de modelo polinomial NARX para o helicóptero
utilizando a resposta em frequência . . . . . . . . . . . . . . . . . . .
94
Comentários Finais . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
97
6 Considerações Finais
6.1
90
99
Sugestões para trabalhos futuros . . . . . . . . . . . . . . . . . . . . . . . . . 101
Referências Bibliográficas
108
Lista de Figuras
3.1
Análise no domínio da frequência de agrupamentos de termos de entrada submetidos a um sinal de frequência única. . . . . . . . . . . . . . . . . . . . . .
31
Análise no domínio da frequência de um modelo com agrupamentos de termos
cruzados submetidos a um sinal de frequência única. . . . . . . . . . . . . . .
34
Análise no domínio da frequência de um modelo com agrupamentos de termos
puros de saída submetidos a um sinal de frequência única. . . . . . . . . . . .
38
4.1
Sistema simulado 1 em representação por blocos. . . . . . . . . . . . . . . . .
50
4.2
Dados gerados pela simulação do sistema 1 para identif cação. . . . . . . . . .
51
4.3
Dados gerados pela simulação do sistema 1 para validação. . . . . . . . . . . .
51
4.4
Ganho das frequências do sistema simulado 1 em dB. . . . . . . . . . . . . . .
52
4.5
Parte dos dados para o sistema 1 excitado a uma frequência de 0,02 rad/s. . . .
53
4.6
DTFT dos dados do sistema simulado 1 excitado a uma frequência de 0,02 rad/s. 54
4.7
Predição livre e análise no domínio da frequência do modelo obtido pelo método
clássico para o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . .
56
Predição livre e análise no domínio da frequência do primeiro modelo obtido
para o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
58
Predição livre e análise no domínio da frequência do segundo modelo obtido
para o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
59
4.10 Predição livre e análise no domínio da frequência do terceiro modelo obtido
para o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
60
3.2
3.3
4.8
4.9
xviii
4.11 Predição livre e análise no domínio da frequência do quarto modelo obtido para
o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
62
4.12 Predição livre e análise no domínio da frequência do quinto modelo obtido para
o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
63
4.13 Evolução do cálculo de ErroW para cada passo para o sistema simulado 1. . .
64
4.14 Parte dos dados gerados pela simulação do sistema 2 para identif cação. . . . .
66
4.15 Parte dos dados gerados pela simulação do sistema 2 para validação. . . . . . .
67
4.16 Ganho das frequências do sistema simulado 2 em dB. . . . . . . . . . . . . . .
67
4.17 Parte dos dados para o sistema 2 excitado a uma frequência de 2,0 × 10−4 rad/s.
69
4.18 Análise no domínio da frequência dos dados do sistema simulado 2 excitado a
uma frequência de 2,0 × 10−4 rad/s. . . . . . . . . . . . . . . . . . . . . . . .
69
4.19 Predição livre e análise no domínio da frequência do modelo obtido pelo método
clássico para o sistema simulado 2. . . . . . . . . . . . . . . . . . . . . . . . .
71
4.20 Evolução do cálculo de ErroW para cada passo para o sistema simulado 2. . .
74
4.21 Predição livre e análise no domínio da frequência do modelo obtido pelo método
de análise da resposta em frequência para o sistema simulado 2. . . . . . . . .
75
4.22 Curva estática do sistema simulado 2 e dos modelos obtidos pelo critério de
ERR e pelo método proposto. . . . . . . . . . . . . . . . . . . . . . . . . . . .
76
5.1
Dados gerados pelo sistema térmico para identif cação. . . . . . . . . . . . . .
81
5.2
Dados gerados pelo sistema térmico para validação. . . . . . . . . . . . . . . .
81
5.3
Ganho das frequências do sistema real térmico em dB. . . . . . . . . . . . . .
82
5.4
Parte dos dados para o sistema térmico excitado a uma frequência de 5,3 × 10−3
rad/s. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
83
Análise no domínio da frequência do sistema térmico excitado a uma frequência
de 5,3 × 10−3 rad/s. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
83
Predição livre e análise no domínio da frequência do modelo obtido pelo método
clássico para o sistema térmico. . . . . . . . . . . . . . . . . . . . . . . . . . .
85
5.5
5.6
xix
5.7
Evolução do cálculo de ErroW para cada passo para o sistema térmico. . . . .
86
5.8
Predição livre do modelo obtido pelo método proposto para o sistema térmico. .
87
5.9
Curva estática do sistema térmico e dos modelos obtidos pelo critério de ERR+AKAIKE
e pelo método proposto. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 88
5.10 Os dois eixos de rotação do helicóptero. . . . . . . . . . . . . . . . . . . . . .
89
5.11 Foto do protótipo de helicóptero com dois GDL. . . . . . . . . . . . . . . . . .
89
5.12 Dados gerados pelo helicóptero para identif cação. . . . . . . . . . . . . . . .
91
5.13 Dados gerados pelo helicóptero para validação. . . . . . . . . . . . . . . . . .
91
5.14 Parte dos dados para o helicóptero excitado a uma frequência de 0,0628 rad/s. .
92
5.15 Análise no domínio da frequência dos dados do sinal de arfagem do helicóptero
excitado a uma frequência de 0,0628 rad/s. . . . . . . . . . . . . . . . . . . . .
92
5.16 Predição livre e análise no domínio da frequência do modelo obtido pelo método
clássico para o helicóptero. . . . . . . . . . . . . . . . . . . . . . . . . . . . .
94
5.17 Evolução do cálculo de ErroW para cada passo para o helicóptero. . . . . . .
95
5.18 Predição livre do modelo obtido pelo método proposto para o helicóptero. . . .
96
Lista de Tabelas
3.1
Resposta em frequência de agrupamentos de termos submetidos a um sinal de
entrada de frequência única ωi . . . . . . . . . . . . . . . . . . . . . . . . . . .
40
4.1
Frequências e respectivos ganhos surgidos na saída do sistema simulado 1. . . .
54
4.2
Frequências e respectivos ganhos surgidos na saída modelo obtido pelo método
clássico para o sistema simulado 1. . . . . . . . . . . . . . . . . . . . . . . . .
56
Frequências e respectivos ganhos surgidos na saída modelo obtido pelo primeiro
passo do método proposto para o sistema simulado 1. . . . . . . . . . . . . . .
58
Frequências e respectivos ganhos surgidos na saída modelo obtido pelo segundo
passo do método proposto para o sistema simulado 1. . . . . . . . . . . . . . .
59
Frequências e respectivos ganhos surgidos na saída modelo obtido pelo terceiro
passo do método proposto para o sistema simulado 1. . . . . . . . . . . . . . .
61
Frequências e respectivos ganhos surgidos na saída modelo obtido pelo quarto
passo do método proposto para o sistema simulado 1. . . . . . . . . . . . . . .
62
Frequências e respectivos ganhos surgidos na saída modelo obtido pelo quinto
passo do método proposto para o sistema simulado 1. . . . . . . . . . . . . . .
63
Frequências e respectivos ganhos surgidos na saída do modelo obtido pelo
método clássico para o sistema simulado 2. . . . . . . . . . . . . . . . . . . .
72
Frequências e respectivos ganhos surgidos na saída do sistema simulado 2. . . .
72
4.10 Valores de ErroW para cada adição de termo no modelo para o sistema simulado 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
73
4.3
4.4
4.5
4.6
4.7
4.8
4.9
xxii
5.1
Frequências e respectivos ganhos surgidos na saída do sistema térmico. . . . .
84
5.2
Frequências e respectivos ganhos surgidos na saída do modelo obtido pelo
método clássico para o sistema térmico. . . . . . . . . . . . . . . . . . . . . .
85
5.3
Valores de ErroW para cada adição e termo no modelo para o sistema térmico.
86
5.4
Frequências e respectivos ganhos surgidos na saída do helicóptero. . . . . . . .
93
5.5
Frequências e respectivos ganhos surgidos na saída do modelo obtido pelo
método clássico para o helicóptero. . . . . . . . . . . . . . . . . . . . . . . . .
94
Valores de ErroW para cada adição e termo no modelo para o sistema térmico.
95
5.6
Lista de Símbolos
F
`
y(k)
u(k)
e(k)
N
τd
ny
nu
ne
c
ψ(k − 1)
T
θ
ˆ
ξ(k)
g
ω
n
nθ
N
ψ
j
x[n]
Ω
A[k]
Gmf
Função
Grau de não linearidade
Sinal de saída
Sinal de entrada
Sinal de ruído
Conjunto dos números naturais
Atraso de transporte
Máximo atraso de saída
Máximo atraso de entrada
Máximo atraso de erro
Parâmetro
Vetor de regressores (pode conter observações até o instante (k − 1))
Transposto
Parâmetros
Valor estimado
Resíduo
Parâmetros estimados
Regressores ortogonais
Número de termos de um modelo
Número de parâmetros de um modelo
Número de amostras
Matriz de regressores
sqrt−1
Sinal
Agrupamento de termos
Peso de senóides
Ganho de uma frequência do modelo
xxiv
Grf
δ
Ganho de uma frequência do sistema
Função impulso unitário
Siglas e Abreviações
AIC
ARX
DTFS
ERR
FOS
FFT
FS
GFRF
MQ
MQE
NARMAX
NARX
OPS
OS
RNA
RMSE
SISO
Critério de Informação de Akaike (Akaike Information Criterion)
Modelo Autoregressivo com Entradas Exógenas (AutoRegressive model with eXogenous inputs)
Séries de Fourier no Tempo Discreto (Discrete Time Fourier Series)
Taxa de Redução de Erro (Error Reduction Ratio)
Busca Ortogonal Rápida (Fast Orthogonal Search)
Transformada Rápida de Fourier (Fast Oourier Transform)
Séries de Fourier (Fourier Series)
Funções de Resposta em Frequência Generalizadas
(Generalized Frequency Response Functions)
Mínimos Quadrados
Mínimos Quadrados Estendidos
Modelo não-linear auto-regressivo de média móvel com entradas exógenas
(Nonlinear AutoRegressive Moving Average model with eXogenous inputs)
Modelo não-linear auto-regressivo com entradas exógenas
(Nonlinear AutoRegressive model with eXogenous inputs)
Busca Ótima de Parâmetros (Optimal Parameter Search)
Busca Ortogonal (Orthogonal Search)
Redes Neurais Artif ciais
Raiz quadrada do erro médio quadrático (Root Mean Square Error)
Uma Entrada e Uma Saída ( Single-Input Single-Onput )
xxvi
Capítulo 1
Introdução
Modelagem matemática de sistemas é uma metodologia muito comum em ciência e na engenharia. Porém mostrava-se um desaf o aos pesquisadores desde as décadas de 50 e 60 pela falta de
ferramentas computacionais para realização de cálculos em série. Mas com o forte avanço tecnológico e o desenvolvimento de recursos computacionais cada vez mais ef cazes, os cálculos
matemáticos já não formariam barreiras para o desenvolvimento da área de modelagem.
Em breves palavras, modelo é uma representação matemática aproximada de um sistema dinâmico real, ou seja, um modelo representa apenas algumas características do sistema real. Nesse
ponto cabe a seguinte pergunta: quais características devem ser reproduzidas pelo modelo? A
resposta a esta pergunta não é única e depende do(s) objetivo(s) para o qual o modelo está sendo
desenvolvido (AGUIRRE, 2004). Uma outra pergunta também poderia ser realizada: Qual(is)
é(são) a(s) região(ões) de interesse em análise deste modelo?
Várias são as representações matemáticas possíveis para sistemas lineares e não lineares. Entre
as representações paramétricas, para um sistema linear existem as funções de transferência, os
modelos em espaço de estados e modelos polinomiais ARX (do inglês - Autoregressive Models with Exogenous Inputs) e ARMAX (do inglês - Autoregressive Moving Average Models
with Exogenous Inputs). Para os sistemas não lineares pode-se citar as representações por Redes Neurais Artif ciais e os modelos polinomiais NARX (do inglês - Nonlinear Autoregressive
Models with Exogenous Inputs) e NARMAX (do inglês - Nonlinear Autoregressive Moving
Average Models with Exogenous Inputs)
De acordo com o nível de conhecimento sobre o processo a ser incorporado ao modelo, os métodos de modelagem podem ser classif cados como em (TULLEKEN, 1993; SJOBERG e outros,
1995):
2
1 Introdução
• modelagem caixa branca (modelo físico ou fenomenológico): nesta abordagem a estrutura e os parâmetros do modelo são determinados por princípios físicos que regem o
comportamento estático e dinâmico do sistema,
• modelagem caixa preta (modelo empírico): a única informação utilizada são dados de entrada e saída do processo. Neste caso, os fenômenos físicos internos ao processo são desconsiderados. Podem ser gerados modelos paramétricos, onde os parâmetros, em geral,
não têm signif cado físico, ou modelos não-paramétricos, representados, por exemplo,
pela resposta ao impulso, resposta ao degrau e resposta em frequência,
• modelagem caixa cinza: nesta técnica, o procedimento de modelagem utiliza tanto dados
de entrada e saída, quanto informações conhecidas do processo. Por exemplo, a estrutura
das equações matemáticas do modelo podem derivar de informações conhecidas a priori
e utiliza-se os dados coletados para identif car os parâmetros deste modelo.
Devido ao conhecimento e tempo necessários para modelar um sistema partindo do equacionamento dos fenômenos envolvidos, nem sempre é viável seguir o procedimento de modelagem
fenomenológica (AGUIRRE, 2004) por isso, em casos onde o processo fenomenológico é deveras complicado, faz-se a opção pela modelagem empírica.
A modelagem empírica é uma alternativa interessante, porque pouco ou nenhum conhecimento
prévio do sistema é necessário na obtenção de um modelo, exigindo-se apenas o conhecimento
dos sinais de entrada e saída do processo. A ciência que estuda modelos matemáticos empíricos
baseados nos sinais de entrada e saída dos sistemas é denominada Teoria de Identificação de
Sistemas.
Uma etapa de fundamental importância na Identif cação de Sistemas é a seleção de estrutura
do modelo. A escolha inadequada pode não incorporar certas características do sistema dinâmico ao modelo. Algumas técnicas são difundidas na literatura, porém, não são ef cientes em
todos os casos e muitas vezes dependem da experiência e da habilidade do desenvolvedor para
mostrarem resultados satisfatórios. Neste trabalho é investigada a seleção de estrutura de modelos matemáticos com o auxílio da resposta em frequência fornecida pelo conjunto de dados
de entrada e saída.
Considerando representações polinomiais, para modelos lineares, a seleção de estrutura limitase a ordem dos polinômios. Em modelos não-lineares é complicado def nir uma ordem signif cativa para a inclusão de termos candidatos no modelo, tendo em vista o elevado número de
termos existentes mesmo para valores relativamente de baixa ordem. Outro problema é def nir
quais serão os cruzamentos de termos visto que tais informações não estão explícitas nos sis-
3
temas ou dados. Aguirre (2004), Aguirre; Billings (1995a) e Mendes; Biillings (2001) já cocomentam sobre as dif culdades na seleção de estruturas de modelos não-lineares.
Um outro problema na identif cação de modelos polinomiais é o fato de que o número de termos candidatos cresce rapidamente com o aumento do grau de não-linearidade e dos máximos
atrasos nos sinais de entrada e saída do modelo. Isso acarreta o aumento da complexidade do
modelo tornando a busca por uma estrutura ideal1 uma tarefa ainda mais árdua.
Alguns métodos para seleção de estrutura foram propostos no decorrer dos estudos sobre a modelagem empírica. Em identif cação caixa-preta, por exemplo, alguns métodos que auxiliam na
seleção de estrutura de modelos não lineares já estão difundidos na literatura dos quais podem
ser citados: critério de informação de Akaike (AKAIKE, 1974; AGUIRRE, 2004), agrupamento
de termos (CORRÊA, 2001; AGUIRRE, 2004), a busca ortogonal, OS, ou busca ortogonal rápida, FOS (do inglês Orthogonal Search e Fast Orthogonal Search) com taxa de redução de erro,
ERR (do inglês Error Reduction Ratio), (KORENBERG, 1985, 1989a,b; BILLINGS; CHEN;
KORENBERG, 1989; CHEN; BILLINGS; LUO, 1989; PAARMANN; KORENBERG, 1992),
a busca ótima de parâmetros com índice de projeção de distância, OPS (do inglês Optimal Parameter Search) (LU; JU; CHON, 2001; ZOU; CHON, 2004) e o SRR (do inglês Simulation
Error Reduction Ratio) (SPINELLI W., 2006; PIRODDI, 2008).
Ressalta-se que a taxa de redução de erro (ERR) ou a busca de parâmetros ótima (OPS), e os
critérios de informação como o critério de Akaike, são complementares entre si. Uma possibilidade, por exemplo, seria o uso da ERR ou OPS para ordenar hierarquicamente os termos
candidatos do modelo. Portanto, tendo-se ordenado os termos, poder-se-ia empregar um critério
de informação para efetuar o corte necessário para a escolha da ordem do modelo representativo. Billings (1980) comenta sobre a dif culdade na busca de representações não lineares e a
impossibilidade de estabelecer apenas uma técnica capaz de fornecer soluções ef cazes em todos
os casos. Zhang (2007) e Wei e Billings (2008) reforçam esta dificuldade.
Um sistema dinâmico pode ser analisado no domínio do tempo e/ou no domínio da frequência.
Por esse motivo, a identif cação deve ser capaz de derivar modelos (lineares ou não-lineares) que
descrevam o comportamento do sistema original no domínio do tempo (equações diferenciais
ou de diferenças) ou no domínio da frequência (resposta em frequência), ponto de interesse
deste trabalho. Seria então possível que a partir de um conjunto de dados de entrada e saída
analisados no domínio da frequência fornecesse informações para a identif cação de modelos
lineares e não-lineares.
1
Uma estrutura ideal para um modelo polinomial é aquela que represente o sistema dinâmico a ser modelado e
tenha a menor complexidade possível.
4
1 Introdução
O termo resposta em frequência signif ca a resposta em regime permanente de um sistema a
uma entrada senoidal. Nos métodos de resposta em frequência, varia-se a frequência do sinal
de entrada dentro de um determinado intervalo e analisa-se a resposta em frequência resultante
(OGATA; KOHN; BARROS MORAIS, 2003). Para sistemas lineares em qualquer uma de suas
faixas de operação, as frequências dos sinais de entrada são as mesmas frequências observadas
nos sinais de saída havendo apenas variações em seus ganhos e/ou fases. O comportamento
dos modelos lineares são como de f ltros que preservam as faixas de frequências de interesse
atenuando aquelas que não explicam o sistema e também inserem ganhos em módulo e avanço
ou atraso de fase nos sinais.
Os sinais de interesse neste trabalho são discretizados, pois são obtidos por placas digitais que
coletam dados em determinados intervalos de tempo. Esse procedimento acarreta o surgimento
de frequências inexistentes no espectro de frequência real dos dados de entrada e saída do sistema gerando o efeito conhecido como aliasing. Neste caso, o que se propõe é que quanto menor
o tempo de amostragem, há maior proximidade do comportamento contínuo. Isso diminuiria
os efeitos da discretização do sinal na análise de resposta em frequência. No entanto, a superamostragem pode provocar outros efeitos. Como exemplo, as altas frequências provocadas por
ruído passam a incorporar o espectro de frequência do sinal causando o mesmo efeito aliasing.
Um outro problema que se manifesta está relacionado à quantidade f nita de dados. As frequências geradas pelo truncamento do sinal podem ser tão expressivas quanto aquelas que realmente
o compõem trazendo dif culdades na análise no domínio da frequência. Por isso, quanto maior
a quantidade de dados, melhor seria a análise no domínio da frequência do sistema. Para isso,
são implementados f ltros conhecidos como janelas que diminuem o efeito do truncamento do
sinal.
A obtenção da resposta em frequência de sinais discretos é amplamente discutida em (HAYKIN;
VAN VEEN, 2001) sendo realizada por exemplo, através de Transformações Discretas de Fourier (FFT do inglês Fast Fourier Transform). No caso das transformadas de Fourier, considerando um sinal como uma superposição ponderada de senóides complexas, aplicando-se
esse sinal a um sistema linear, a saída do sistema será uma superposição ponderada das respostas do sistema a cada senóide complexa. Porém, em sistemas não lineares o mesmo não
ocorre, podendo haver surgimento e/ou desaparecimento nas saídas de frequências aplicadas
na excitação da entrada. Os sistemas não-lineares são conhecidos pela presença de harmônicos, inter-modulações complexas e até mesmo o caos, onde há transferência de energia entre
as diferentes frequências fornecendo resultados no domínio da frequência bem diferentes das
frequências de entrada (RUGH, 1981; LANG e outros, 2007; WU; LANG; BILLINGS, 2007;
JING; LANG; BILLINGS, 2008).
5
As respostas em frequências de modelos não lineares podem ser obtidas através das Funções
de Resposta em Frequência Generalizadas (GFRF - do inglês Generalized Frequency Response
Functions) que estendem o conceito de funções de resposta em frequência lineares a uma ampla
classe de sistemas não-lineares que possuem memória e podem ser descritos pelo modelo de
séries de Volterra (VOLTERRA, 1930; BOYD; CHUA; DESOER, 1984).
No entanto, o que se deseja é obter o modelo matemático a partir do conjunto de dados de
entrada e saída do sistema e com pouco ou nenhum outro conhecimento prévio. O que pode
ser feito com esses dados são as suas análises no domínio da frequência através da FFT. Tal
ferramenta permite visualizar qual foi o conjunto de frequências existente nos dados de entrada,
ou seja, aquelas frequências às quais o sistema foi submetido, e as frequências que compõem
o sinal gerado em sua saída. Isso é de fundamental importância pois a partir dessa análise já é
possível identif car possíveis não linearidades do sistema modelado. A partir de então pode ser
considerado que a análise da frequência da saída gerada por um determinado sistema não-linear,
além das informações dos dados de entrada, possuem informações sobre o sistema em que se
aplica tal entrada.
O método de análise de resposta em frequência de sistemas não-lineares baseados na teoria de
séries de Volterra foi iniciada em na década de 1950 quando o conceito de GFRF de sistemas não
lineares foi apresentado. As GFRF’s foram def ninas como transformações multidimensionais
de Fourier dos kernels (do inglês Núcleos) de Volterra nas expansões de sistemas não-lineares
nas séries de Volterra (LANG; BILLINGS, 2000).
Lang e outros (2007) realizaram um trabalho que procurava investigar as frequências de saída
de sistemas não-lineares com o objetivo de contornar algumas dif culdades na análise e desenvolvimento no domínio da frequência. Ainda propuseram uma metodologia para obter as
respostas em frequência de modelos polinomiais de equações de diferenças a partir dos termos
do modelo e de seus coeficientes. Wu; Lang; Billings (2007) comentam sobre o desenvolvimento de um algoritmo capaz de determinar as faixas de frequência de saída de sistemas nãolineares.
No entanto, poucos estudos baseados no domínio da frequência focaram-se na seleção de estruturas de modelos polinomiais discretos. A grande maioria envolveu a determinação das GFRF’s
a partir de modelos já obtidos e parametrizados.
6
1 Introdução
1.1 Motivação
A motivação para realização deste trabalho está no fato de que o processo de escolha de estrutura para modelos polinomiais ainda encontra-se em aberto. Pouco tem sido discutido quanto
a relação de análise de resposta em frequência de um sistema e os termos necessários para
reproduzi-la utilizando um modelo polinomial discreto. Logo, neste trabalho são investigadas
as possibilidades de seleção de termos e consequentemente seleção de estruturas, a partir análise
no domínio da frequência do conjunto de dados de entrada e saída do sistema.
1.2 Objetivos
Os objetivos deste trabalho são:
1. Investigar a relação entre a inclusão de combinações dos vários termos candidatos de um
modelo discreto de equações de diferenças NARX e as mudanças espectrais dos sinais de
saída em relação às mudanças espectrais dos sinais de entrada;
2. Propor um procedimento de auxílio na escolha da estrutura (grau de não linearidade e
termos candidatos) de um modelo de equações de diferenças que se baseie nas análises
do domínio da frequência do conjunto de dados de entrada e saída;
3. Comprovar a ef ciência do método proposto através de testes em dois sistemas simulados
e em dois sistemas reais.
1.3 Organização do texto
Esta dissertação está organizada em seis capítulos da seguinte forma:
Capítulo 2: Fundamentação teórica. Neste capítulo, a fundamentação teórica para a realização deste trabalho é apresentada. A Identif cação de Sistemas é apresentada em cada uma
de suas etapas. Em seguida são mostradas as técnicas de análise de resposta em frequência
utilizadas na aplicação do método proposto.
Capítulo 3: Auxílio da seleção de estruturas a partir das respostas em frequências dos
dados de entrada e saída. Neste capítulo, procedimentos de auxílio na seleção de estrutura de
modelos de equações de diferenças a partir das respostas em frequência dos dados de entrada e
1.3 Organização do texto
7
saída são propostos e discutidos.
Capítulo 4: Estudo em sistemas simulados. Neste capítulo, os procedimentos apresentados
no capítulo 4 são aplicados à dois sistemas simulados. O primeiro sistema simulado apresenta
uma não linearidade mais branda. O segundo sistema simulado apresenta um grau de não linearidade mais severa. São apresentados e discutidos os resultados para estes dois casos.
Capítulo 5: Estudo em sistemas reais. Este capítulo apresenta a aplicação do procedimento
proposto em dois sistemas reais. O primeiro composto de um forno de bancada. E o segundo
um helicóptero didático de dois graus de liberdade. São apresentados e discutidos os resultados
para os dois casos.
Capítulo 6: Considerações finais. Este capítulo faz as considerações f nais a respeito do
trabalho realizado e apresenta sugestões para trabalhos futuros.
8
1 Introdução
Capítulo 2
Fundamentação Teórica
2.1 Introdução
Neste capítulo será apresentada uma breve fundamentação teórica que sustenta os estudos realizados neste trabalho. Num primeiro momento são apresentadas as técnicas de identif cação,
suas etapas e métodos mais utilizados. Depois são apresentadas as técnicas de análise no domínio da frequência.
2.2 Identificação de sistemas
O comportamento de um sistema pode ser representado por diferentes tipos de modelos. A escolha de um modelo para um determinado sistema depende tanto do foco de interesse de estudo
quanto da aplicação deste modelo. Por exemplo, compreender certas dinâmicas do sistema sob
diversas condições de operação, por exemplo, em um sistema industrial dif cilmente seria possível excitar uma planta em toda sua faixa de operação pois poderiam causar danos ao sistema
ou até mesmo perdas materiais. Uma outra necessidade seria predizer o comportamento de um
determinado sistema, tal como a prevenção de catástrofes naturais evitando perda de vidas. Outras aplicações seriam: otimizar o comportamento do sistema, analisar e projetar controladores,
permitir detecção ef ciente de falhas no sistema, estimar variáveis do processo que não podem
ser medidas diretamente, permitir o estudo do sistema em regiões de operação despendiosas
ou problemáticas no sistema real, permitindo um treinamento de operação seguro e ef ciente
(MATKO; ZUPANCIC; KARBA, 1992).
10
2 Fundamentação Teórica
O comportamento de um sistema pode ser representado por diferentes tipos de modelos. A
escolha do tipo de modelo é inerente a característica a qual deseja-se que o modelo represente
assim como o propósito ou o objetivo para o qual o modelo será desenvolvido. Gevers (2005)
afirma que se o modelo é somente uma aproximação do sistema real, então a qualidade do
modelo deve ser dependente da aplicação desejada.
Dentre esses modelos, os matemáticos permitem aplicações avançadas, sendo, por esse motivo, os mais utilizados, seja na engenharia, biologia, medicina, economia ou em outras áreas
(LJUNG, 1999).
A modelagem matemática constitui-se de duas vertentes principais. A primeira delas é a modelagem pela física do processo também conhecida como modelagem caixa branca. Modelagem
esta que procura, através de relações físico-matemáticas, explicar o comportamento de um determinado sistema. Muitas vezes torna-se deveras complexa pois nem sempre em todo o comportamento do sistema há relações físicas bem def nidas ou tempo suf ciente para desenvolver
um modelo a partir de todas as equações que regem a física do processo. Devido a esta demanda criada pela impossibilidade de modelar fenomenologicamente todos os sistemas, surgiu
outra vertente: A Identificação de Sistemas, ciência que procura explicar os fenômenos de um
processo através apenas de observações, ou seja, a partir dos dados de entrada e saída por isso
conhecida também como modelagem caixa-preta.
A identif cação de sistemas apresenta uma grande vantagem pela relativa facilidade de obtenção
de modelos em virtude do pouco ou nenhum conhecimento prévio do sistema exigido. A
desvantagem até o momento se restringe ao número excessivo de parâmetros e a ausência de
signif cado físico dos modelos obtidos.
Neste capítulo é apresentada uma breve discussão referente as principais etapas da identif cação
de sistemas.
2.2.1 Etapas da Identificação
As principais etapas da Identif cação de Sistemas são:
• Testes dinâmicos e a coleta de dados;
• escolha de representação matemática a ser usada;
• determinação da estrutura do modelo;
• estimação dos parâmetros;
2.2 Identif cação de sistemas
11
• validação do modelo.
Nas próximas seções as etapas supracitadas serão abordadas de forma sucinta. Outros detalhes,
aplicações, exemplos e discussões podem ser encontrados em (AGUIRRE, 2004).
2.2.1.1 Testes dinâmicos e a coleta de dados
Três aspectos precisam ser analisados: onde excitar o sistema, que tipos de sinais usar a f m de
obter dados que sejam representativos dessa dinâmica e como amostrar tais dados.
O sinal de excitação apresenta duas importantes funções. A primeira delas é conseguir excitar todas as características dinâmicas e estáticas do sistema em toda a faixa de frequência de
interesse. A segunda está relacionada ao perf l de amplitudes desse sinal, que é responsável
pela excitação das não-linearidades presentes no sistema. Características dinâmicas e estáticas
que não forem excitadas não aparecerão nos modelos. O que não estiver nos dados não será
identif cado.
Na identif cação de sistemas não-lineares é comum usar um gerador de números aleatórios para
gerar sinais de entrada persistentemente excitantes de ordem elevada. Tendo def nido o tempo
de amostragem, mantém-se constante cada valor escolhido aleatoriamente por um tempo em
torno de 3 a 5 intervalos de amostragem.
Uma etapa importante no processo de identif cação é a escolha do tempo ou período de amostragem, que corresponde ao intervalo entre duas amostras. Na prática, a frequência de amostragem
é normalmente escolhida entre 5 e 10 vezes a maior frequência de interesse.
Dois problemas podem surgir no mal dimensionamento do tempo de amostragem. O primeiro
deles refere-se a subamostragem, onde o tempo de amostragem é muito grande e causa mal
representação da dinâmica real do sistema. O outro problema é o de super-amostragem, onde
o tempo de amostragem é muito pequeno, que como descrito por (BILLINGS; AGUIRRE,
1995), um sinal super-amostrado, provoca instabilidade numérica e elevado esforço computacional devido ao mal condicionamento da matriz de regressores, dif cultando a determinação da
estrutura do modelo. Quando um sinal é super-amostrado, pode-se decimá-lo a f m de que se
torne devidamente amostrado.
2.2.1.2 Escolha da representação matemática
Grandes são as possibilidades de representações não-lineares na identif cação de sistemas. A
escolha do tipo de representação depende de fatores já discutidos na seção 2.1 dos quais pode ser
acrescentada a disponibilidades das ferramentas de implementação da representação escolhida.
12
2 Fundamentação Teórica
Porém na prática as representações são selecionadas de acordo com a conveniência do modelador, a complexidade da aplicação ou a possibilidade de interligação com outros conceitos.
As representações não-lineares que mais se destacam, podem ser citadas, as redes neurais artif ciais (RNA), defendidas e melhor discutidas em (BRAGA; CARVALHO; LUDEMIR, 2000;
AMARAL, 2001; HAYKIN, 2001). Tais representações são modelos matemáticos inspirados
no funcionamento do cérebro humano, formado por um conjunto de células chamadas neurônios
fazendo um paralelo ao grupo de unidades (neurônios) que em conjunto, formam um sistema
inteligente capaz de representar um sistema (cérebro).
Uma outra representação matemática que pode ser citada é a série de Volterra (VOLTERRA,
1930; BILLINGS, 1980). A utilização de uma série de Volterra permite caracterizar um processo não-linear de forma a se ter uma idéia do seu signif cado físico (RUGH, 1981).
Entre as representações não-lineares, pode-se destacar também os modelos de blocos interconectados (WIENER, 1958; WIGREN, 1993; PEARSON; POTTMANN, 2000; COELHO,
2002), compostos por um modelo dinâmico linear em cascata com uma função estática nãolinear. Quando a função estática precede o modelo dinâmico, tem-se o modelo de Hammerstein. Por outro lado, quando o modelo dinâmico precede a função estática, tem-se o modelo de
Wiener.
Um outro tipo de representação não-linear são os modelos polinomiais e racionais (BILLINGS;
CHEN, 1989a,b; BILLINGS; TAO, 1991; BILLINGS; ZHU, 1993, 1994; CORRÊA, 2001;
CAMPOS, 2007).
Na área de identif cação de modelos não lineares, o modelo NARMAX polinomial (LEONTARITIS; BILLINGS, 1985a,b). Tais modelos são extensões dos modelos NARX porém com
a inclusão de termos de ruído para evitar polarização1 de parâmetros.
O modelo NARMAX polinomial é um modelo discreto que explica o valor de saída em função
de valores prévios dos sinais de saída y(k), de entrada u(k) e de ruído e(k) e é representado da
seguinte forma:
y(k) = F ` y(k − 1), . . . ,y(k − ny ),u(k − τd ), . . . ,u(k − τd − nu + 1),
e(k − 1), . . . ,e(k − ne ) + e(k),
1
(2.1)
A polarização é o desvio entre o valor esperado de uma determinada variável aleatória (a variável estimada) e
uma variável determinística (o valor real), ainda que não seja conhecida.
2.2 Identif cação de sistemas
13
sendo que e(k) indica os efeitos que não podem ser bem representados por F ` [·], que é uma
função polinomial com grau de não-linearidade ` ∈ N. Os termos τd ,ny ,nu e ne representam,
respectivamente, o atraso de transporte e os máximos atrasos em y, em u e em e.
A parte livre do ruído da equação (2.1) pode ser expandida como o somatório de termos com
graus de não-linearidade variando na faixa 1 ≤ m ≤ `. Assim, cada termo de grau m poderá
conter um fator de grau p do tipo y(k − i) e um fator de grau (m − p) do tipo u(k − i) sendo
multiplicado por um parâmetro representado por cp,m−p (n1 , . . . ,nm ). Matematicamente, tem-se
(PEYTON-JONES; BILLINGS., 1989):
y(k) =
y ,nu
m nX
X̀ X
cp,m−p (n1 , . . . ,nm )
m=0 p=0 n1 ,nm
p
Y
y(k − ni )
i=1
m
Y
u(k − ni ),
(2.2)
i=p+1
sendo que
ny ,nu
X
n1 ,nm
≡
ny
X
n1 =1
···
nu
X
.
(2.3)
nm =1
O limite superior será ny se o somatório se referir a fatores do tipo y(k − ni ). Por outro lado,
se o somatório se referir a fatores do tipo u(k − ni ), tal limite será nu .
2.2.1.3 Seleção da Estrutura do Modelo
Parte crucial na Identif cação de Sistemas, a seleção de estrutura ainda é um problema. Várias
técnicas foram propostas a f m de selecionar os melhores regressores de modelos NARX/NARMAX embora ainda não se tenha uma solução que def na com ef ciência a melhor estrutura para
cada caso identif cado.
Alguns casos de seleção de estrutura que obtiveram melhores desempenhos na Identif cação
de Sistemas não lineares foram propostos em (KORENBERG e outros, 1988; LEONTARITS; BILLINGS, 1987; BILLINGS; CHEN; KORENBERG, 1989; AGUIRRE; BILLINGS,
1995a; MAO; BILLINGS, 1997);(MENDES; BILLINGS, 2001; PIRODDI; SPINELLI, 2003;
RODRÍGUEZ-VÁZQUEZ e outros, 2004; ZHANG; ZHU; LONGDEN, 2007; WEI; BILLINGS,
2008).
Para modelos polinomiais a seleção de estrutura trata da def nição do grau de não-linearidade
` e dos máximos atrasos em y, em u e em e. É importante salientar que quanto maior o grau
de não linearidade e os máximos atrasos, maior será a variedade de termos possíveis multipli-
14
2 Fundamentação Teórica
cando a possibilidade de combinações e consequentemente de estruturas possíveis. Aguirre;
Billings (1995b) chamam atenção para o fato que se o número de termos for maior que o
número de termos necessários, ou seja, um modelo sobreparametrizado, pode acarretar instabilidade numérica na estimação dos parâmetros (seção 2.2.1.6), esforço computacional além do
necessário e fazer com que o modelo represente dinâmicas não existentes no sistema representado.
Um outro cuidado que o modelador deve ter é que na tentativa de melhor explicar os dados
de identif cação do sistema, os modelos tendem a se tornar excessivamente complexos. A utilização desses modelos na prática torna-se inviável, pois muito tempo e esforço computacional
serão necessários para implementar a estimação de parâmetros, simular o modelo e predizer o
comportamento do sistema.
Neste trabalho considerou-se que o critério ERR (do inglês error reduction ratio) foi usado para
ordenar os termos candidatos do modelo por ordem de importância na explicação da dinâmica
do sistema e tendo-se ordenado os termos, empregou-se o critério de informação de Akaike
para efetuar o corte necessário para a escolha da ordem do modelo representativo para comparação com a técnica investigada. Tanto o critério de Akaike quanto o ERR serão brevemente
explanados nas seções 2.2.1.5 e 2.2.1.4
2.2.1.4 Critério ERR
O critério ERR baseia-se na redução do erro de predição, avalia a importância dos termos do
modelo em função da sua habilidade de explicar a variância do sinal de saída.
A f m de def nir a taxa de redução de erro (ERR), será def nido o seguinte modelo NARMAX
geral:
y(k) = ψ T (k − 1)θ̂ + ξ(k)
nθ
X
θ̂i ψi (k − 1) + ξ(k),
=
i=1
sendo que o termo nθ corresponde ao número de parâmetros.
Analogamente, um modelo representado em base ortogonal pode ser expresso por:
(2.4)
2.2 Identif cação de sistemas
15
y(k) =
nθ
X
ĝi wi (k − 1) + ξ(k),
(2.5)
i=1
sendo que ĝi representa os parâmetros estimados e wi , os regressores ortogonais sobre os dados.
Na ausência de regressores, ou seja, quando nθ = 0, o erro de predição corresponde ao próprio
sinal de saída y(k).
A taxa de redução de erro devido à inclusão do i-ésimo regressor no modelo é def nida como
(CHEN; BILLINGS; LUO, 1989):
[ERR]i =
ĝi2 hwi ,wi i
,
hy,yi
(2.6)
sendo que o ERR indica a porção da variância da saída explicada pela inclusão de um novo
termo no modelo. Normalmente, a escolha dos termos candidatos ao modelo é feita considerando os regressores que apresentam os maiores valores de ERR.
2.2.1.5 Critério de Akaike
O método mais utilizado para estimar o número de termos que devem ser incluídos em modelos
dinâmicos é o critério de informação de AKAIKE (1974). De acordo com o critério de Akaike,
o número ótimo de termos em um modelo deve minimizar a seguinte função de custo:
AIC = N ln(V ar[ξ(k)]) + 2 nθ ,
(2.7)
sendo que N corresponde ao número de amostras, V ar[ξ(k)], à variância do resíduo e nθ , o
número de parâmetros.
A função representada na equação (2.7) signif ca que à medida que os termos são incluídos no
modelo, o número de graus de liberdade aumenta permitindo um ajuste dos dados mais exato.
Assim, V ar[ξ(k)] diminui à medida que nθ aumenta. A partir de um determinado momento, a
redução na variância do resíduo devido à inclusão de um novo termo passa a ser insignif cante
e, consequentemente, não justif ca mais a inclusão deste último termo. Neste ponto é efetuado
o corte.
16
2 Fundamentação Teórica
2.2.1.6 Estimação de Parâmetros
Determinada a estrutura do modelo, deve-se estimar os parâmetros correspondentes a cada um
de seus termos. Para isso, utiliza-se do conjunto de dados obtidos ou simulados separados para
esta fase do procedimento. Tais dados são comumente chamados de dados de identif cação.
Uma outra parcela dos dados é separada para a realização da validação do modelo que será
melhor discutida na seção 2.2.1.7.
Representando a equação (2.1) na forma:
y(k) = ψ T (k − 1)θ̂ + ξ(k)
= ŷ(k) + ξ(k),
(2.8)
sendo que k indica o instante considerado, ψ(k − 1) corresponde ao vetor de regressores, que
pode conter observações até o instante (k − 1), θ̂ representa o vetor de parâmetros estimados2 e
ξ(k), o erro cometido pelo modelo ao tentar explicar y(k) como ψ T (k − 1)θ̂.
Tomando a equação (2.8) ao longo do conjunto de dados e representando o conjunto de equações
resultante em forma matricial, tem-se:
y = Ψθ̂ + ξ,
(2.9)
sendo que Ψ representa a matriz de regressores e ξ representa o erro cometido pelo modelo ao
tentar explicar y(k).
Isolando-se ξ na equação (2.9) tem-se:
ξ = y − Ψθ̂.
(2.10)
No estimador de mínimos quadrados, os parâmetros são estimados minimizando-se a seguinte
função de custo:
JMQ =
N
X
ξ(i)2 = ξ T ξ.
i=1
2
O símbolo ( ˆ ) sobre as variáveis faz referência a valores estimados.
(2.11)
2.2 Identif cação de sistemas
17
Substituindo a equação (2.10) em (2.11), tem-se:
θ̂ MQ = (ΨT Ψ)−1 ΨT y.
(2.12)
que é a expressão do estimador de mínimos quadrados (MQ).
Em algumas situações, o MQ torna-se polarizado, levando ao uso de soluções como o estimador
de mínimos quadrados estendido (MQE).
O MQE consiste em estender a matriz de regressores Ψ incluindo regressores de ruído para que
a parte modelável do erro seja modelada tornando o ruído ξ branco. Este procedimento faz com
que não haja correlação entre a matriz estendida de regressores e o resíduo. Para isso, o resíduo
é utilizado como estimativa do ruído usado na formação da nova matriz de regressores.
2.2.1.7 Validação do modelo
Após a obtenção do modelo utilizando-se das quatro etapas anteriores, é necessário verif car se
este atende às expectativas do modelador. O modelo que atenderá as expectativas será aquele
que suprir a necessidade para a qual foi projetado e será usado. Como já foi dito, é impossível
que o modelo represente todas as dinâmicas do sistema. Porém um modelo será considerado
válido se incorporar as características do sistema que são fundamentais para a aplicação desejada.
Aguirre (2004) destaca que um cuidado básico é não usar os dados utilizados para obter o
modelo na validação. A motivação para tal cuidado é que a validação procura saber a capacidade
de generalização do modelo e não somente o ajuste aos dados usados na identif cação. Essa
capacidade de generalização signif ca que o modelo realmente incorporou a dinâmica do sistema
e não apenas ajustou aos dados.
Uma técnica muito utilizada para validação de modelos é a predição livre, conhecida também
como predição de inf nitos passos à frente, no qual o modelo é simulado indef nidamente retroalimentado com predições passadas. Porém, para iniciar o processo, são necessários dados
medidos retirados da parcela separada para validação.
Para quantif car a qualidade de modelo e obter um parâmetro de comparação de desempenho o
RMSE (do inglês Root Mean Squared Error) é largamente utilizado e é dado por:
18
2 Fundamentação Teórica
RMSE =
qP
2
y(k) − ŷ(k)
qP
2 ,
N
y(k)
−
y
k=1
N
k=1
(2.13)
sendo que ŷ(k) é a predição livre do sinal e y é o valor médio do sinal medido y(k). É importante
ressaltar que quanto menor for o valor do RMSE, maior será a capacidade do modelo de se
ajustar aos dados. Se forem utilizados os dados de validação como exige a técnica, comprovase que quão menor seja o índice RMSE, maior será a capacidade de generalização do modelo.
Até então, tratou-se da validação dinâmica, ou seja utilizando pares de dados de entrada e saída
transcorrentes no tempo. Porém, existe também a validação estática, ou seja, aquela que usa
a relação entre os dados de entrada e saída sem a inf uência do estado de transição temporal
também conhecido como estado estacionário. Corrêa (2001), Coelho (2002) e Campos (2007)
utilizaram desta técnicas em seus trabalhos.
2.3 Transformada de Fourier
A quase totalidade dos fenômenos físicos de natureza periódica tem sua origem em ondas (sonoras, elétricas, hidráulicas e mecânicas) que apresentam formas bem def nidas.
São quatro as representações de Fourier. Cada uma aplicável a uma classe diferente de sinais.
Essas quatro classes são def nidas pelas propriedades de periodicidade de um sinal mostrando
se ele é contínuo ou discreto (HAYKIN; VAN VEEN, 2001).
Os sinais periódicos podem ser representados como séries de Fourier FS (do inglês Fourier Series). Essas séries aplicam-se a sinais periódicos no tempo contínuo. Também existem as séries
de Fourier de tempo discreto (DTFS do inglês Discrete Time Fourier Series) que se aplicam-se
a sinais periódicos discretizados. Sinais não-periódicos contínuos são representados pela Transformada de Fourier, FT (do inglês Fourier Transform). E os sinais não-periódicos discretizados
são representados pela Transformada de Fourier de Tempo Discreto, DTFT (do inglês Discrete
Time Fourier Transform).
Considerando que os sinais utilizados neste trabalho são sinais discretizados e os sinais de
interesse são sinais não-periódicos, a discussão nesta seção se limitará às DTFT’s.
Na seção 2.3.1 será apresentada uma avaliação da base matemática para converter qualquer sinal
não-periódico nas suas componentes de frequência utilizando a sua DTFT. Posteriormente, na
seção 2.3.2 serão apresentados os efeitos negativos da amostragem dos sinais e o que pode ser
2.3 Transformada de Fourier
19
feito para amenizá-los.
2.3.1 As Transformadas de Fourier no tempo discreto
Nesta discussão, empregar-se-á uma abordagem intuitiva procurando simplif car a demonstração
da DTFT. Partindo-se da DTFS, descreve-se um sinal não-periódico como o limite de um sinal
periódico cujo período, N, aproxima-se do inf nito. Supõe-se que o sinal não-periódico seja
representado por um período simples do sinal periódico que está centralizado na origem, e o
que o limite, quando N se aproxima do inf nito, é tomado de maneira simétrica. Admitindo
que x̃[n] seja um sinal periódico com período N = 2M + 1. Def ne-se, portanto, o sinal não
periódico de duração f nita x[n] como um período de x̃[n], como é mostrado por
x[n] =
x̃[n], −M ≤ n ≤ M
0, |n| > M
(2.14)
Procurando replicar o sinal não-periódico inf nitas vezes para longe da origem tanto pela direita
como pela esquerda tem-se M → ∞ e desta forma pode-se expressar
x[n] = lim x̃[n].
M →∞
(2.15)
Partindo da representação de DTFS do sinal periódico x̃[n] tem-se
x̃[n] =
M
X
X[k]ejkΩ0 n
(2.16)
k=−M
sendo Ω0 = 2π/N a frequência fundamental de x[n] que é o vetor de dados. A frequência da
k-ésima senóide na superposição é kΩ0 . X[k] é o peso aplicado à k-ésima senóide complexa e
é dado por
M
X
1
x̃[n]e−jkΩ0 n .
X[k] =
2M + 1 n=−M
(2.17)
Uma vez que x̃[n] = x[n] para −M ≤ n ≤ M, a equação (2.17) pode ser reescrita da seguinte
forma:
20
2 Fundamentação Teórica
∞
X
1
X[k] =
x[n]e−jkΩ0 n ,
2M + 1 n=−∞
(2.18)
partindo do fato que x[n] = 0 para |n| > M. Def ne-se agora uma função contínua de frequência X(ejΩ ) cujas amostras em kΩ0 são iguais aos coef cientes da DTFS normalizados por
2M + 1. Logo,
∞
X
X(ejΩ ) =
x[n]e−jΩn ,
(2.19)
n=−∞
de forma que X[k] = X(ejkΩ0 )/(2M + 1). Substituindo esta def nição de X[k] na equação
(2.16) tem-se:
M
X
1
x̃[n] =
X(ejkΩ0 )ejkΩ0n
2M + 1 k=−M
(2.20)
e usando a relação Ω0 = 2π/(2M + 1), escreve-se:
M
1 X
x̃[n] =
X(ejkΩ0 )ejkΩ0n Ω0 .
2π k=−M
(2.21)
Considerando x[n] como o valor limite de x̃[n] quando M → ∞, avaliar-se-á o efeito de
M → ∞ sobre a frequência fundamental Ω0 . À medida que M se eleva, Ω0 decresce e os
espaçamentos entre os harmônicos da DTFS diminui. Isso pode ser notado matematicamente
através da equação (2.19). Combinando a equação (2.15) com a equação (2.21) tem-se:
M
1 X
x[n] = lim
X(ejkΩ0 )ejkΩ0 n Ω0 .
M →∞ 2π
k=−M
(2.22)
Aproximando a equação (2.23) pela regra retangular para uma integral e lembrando que Ω =
kΩ0 , de forma que dΩ = Ω0 , tem-se:
1
x[n] =
2π
Z
π
−π
X(ejΩ )ejΩn dΩ.
(2.23)
2.3 Transformada de Fourier
21
Dessa forma x[n] é expresso como uma superposição ponderada de senóides de tempo discreto.
Neste caso, a superposição é uma integral e a ponderação em cada senóide é 12 πX(ejΩ)dΩ.
Logo, a representação por DTFT é expressa como:




x[n] =
1 Rπ
X(ejΩ )ejΩn dΩ
2π −π
(2.24)


 X(ejΩ ) = P∞
−jΩn
n=−∞ x[n]e
Logo, pode-se dizer que X(ejΩ ) e x[n] são um par de DTFT, descrito por:
FT
−−−−−→ X(ejΩ ),
x[n]←−−DT
ou seja, a transformada X(ejΩ ) descreve o sinal x[n] como uma função de frequência senoidal
Ω e é denominada representação no domínio da frequência de x[n]. A amostragem de um sinal
no tempo gera no domínio da frequência réplicas deslocadas do espectro do sinal original.
Mais informações sobre a Transformada de Fourier de Tempo Discreto pode ser encontrada em
(HAYKIN; VAN VEEN, 2001).
2.3.2 Efeitos da amostragem de sinais
2.3.2.1 Aliasing
A operação de amostragem gera um sinal de tempo discreto a partir de um sinal de um tempo
contínuo, ou seja, representa-se matematicamente o sinal amostrado como o produto de um sinal
contínuo original por um trem de impulsos. Os efeitos da subamostragem e super-amostragem
no processo de identif cação foram brevemente discutidos na seção 2.2.1.1. No entanto, quando
um sinal é amostrado, é necessário que possa se reconstruir o sinal original (contínuo) a partir
do sinal digitalizado. Isso garante que todas as características do sinal original inclusive sua
Transformada de Fourier sejam preservadas.
Quando há subamostragem do sinal em estudo, pode-se haver um efeito chamado aliasing,
que em português signif ca disfarce. Esse disfarce ocorre quando a DTFT deixa de ter uma
correspondência biunívoca com o sinal de tempo contínuo original. Isso signif ca que não é
possível usar o espectro do sinal amostrado para analisar o sinal no tempo contínuo, e que não
se pode reconstruir de forma única o sinal de tempo contínuo a partir de suas amostras.
22
2 Fundamentação Teórica
O aliasing distorce o espectro do sinal original. Há uma superposição entre o espectro original e
suas réplicas acarretando o fenômeno de um componente de alta frequência assumir a identidade
de um de baixa frequência. Isto sugere que uma condição para a correspondência única entre
o sinal de tempo contínuo e suas amostras é equivalente à condição para o impedimento do
aliasing. Esta exigência é formalmente estabelecida da seguinte maneira:
Teorema 2.3.1 (Teorema da Amostragem) Admitindo-se que x(t)←−−F−T−−→ X(jω) represente um
sinal de faixa limitada, de forma que X(jω) = 0 para |ω| > ωm , sendo ω uma frequência qualquer do sinal e ωm a mais alta frequência do sinal. Se ωs > 2ωm , em que ωs = 2π/τ é a
frequência de amostragem e τ , o tempo de amostragem, então x(t) é determinado de maneira
única por suas amostras x(nτ ), com n = 0, ± 1, ± 2, . . . .
A frequência mínima de amostragem, 2ωm , é denominada taxa de amostragem de Nyquist ou
taxa de Nyquist.
Isso leva a um relaxamento em escolher uma taxa de amostragem superior à taxa de Nyquist.
Porém, embora o sinal no qual se tem interesse possa ter uma frequência máxima ωm , em geral,
o sinal de tempo contínuo terá energia em frequências mais elevadas devido à presenta de ruído
ou de outras características não ideais. Portanto, uma super-amostragem também pode provocar
aliasing na DTFT de um sinal de interesse.
2.3.2.2 Janelamento do sinal
Em aplicações práticas envolvendo a amostragem de sinais pode-se obter somente uma gravação
f nita do sinal. Isso resulta em uma forma de onda truncada que possui características espectrais diferentes do sinal original. Tal descontinuidade produz a perda da informação espectral
original.
Uma maneira de aumentar as características espectrais de um sinal amostrado é a aplicação
de janelas sobre o mesmo. Ao analisar uma sequência de dados f nita através da DTFT, o
janelamento minimiza as margens de transição em formas de onda truncadas, reduzindo dessa
forma a perda espectral.
O janelamento de um sinal no domínio do tempo é equivalente um f ltro passa baixas. Resumidamente, o que se faz é multiplicar o sinal f nito pela função que representa a janela. Devido
a multiplicação no domínio do tempo ser equivalente à convolução no domínio da frequência,
o espectro de um sinal ajanelado é a convolução do espectro do sinal original com o espectro
da janela. Dessa maneira, o janelamento modif ca a forma do sinal tanto no domínio do tempo
quanto no da frequência.
2.4 Considerações Finais
23
Os principais tipos de janelas disponíveis são citados a seguir:
• Retangular: Aplicar uma janela retangular é equivalente a não utilizar qualquer janela.
A janela retangular possui o maior volume de perda espectral. Possui o valor igual a 1
sobre todo o seu intervalo de tempo.
• Hanning: Esta janela possui uma forma similar aquela de meio ciclo de uma forma de
onda cossenoidal. Uma janela de tamanho N está def nida através da Equação :
w[n] = 0,5 − 0,5cos(2πn/N), para n = 0,1,2, . . . , N − 1
(2.25)
• Hamming: Essa janela é uma versão modif cada da janela de Hanning. Sua forma também é similar a de uma onda cossenoidal. Uma janela de tamanho N é def nida pela
Equação .
w[n] = 0,54 − 0,46cos(2πn/N), para n = 0,1,2, . . . , N − 1
(2.26)
Outras janelas também podem ser citadas tais como a janela de Kaiser-Bessel, Triangular, Flattop, Exponencial, entre outras. A escolha correta da janela depende de um certo conhecimento
a priori do sinal em análise. Como em grande parte das vezes isso não é possível, deve-se experimentar os diversos métodos de janelamento existentes e encontrar o que melhor se adapta a
análise em questão.
2.4 Considerações Finais
Neste capítulo foram discutidas as principais etapas da identif cação de sistemas. Métodos de
tratamento de dados, seleção de estruturas e estimação de parâmetros foram apresentados e
discutidos. Alguns trabalhos realizados utilizando tais técnicas foram citados. Mais detalhes
sobre identif cação de sistemas podem ser vistos em (AGUIRRE, 2004).
Também foram apresentadas a def nição e a forma de obtenção da Transformada de Fourier de
sinais discretos. Alguns problemas que podem deturpar a transforma também foram apresentados e brevemente discutidos. Mais detalhes sobre essa discussão, conceitos e aplicações podem
ser encontrados em (HAYKIN; VAN VEEN, 2001).
Os métodos descritos nesse capítulo foram utilizados para identif cação de sistemas de dados
usados neste trabalho e para obtenção da resposta em frequência de sinais aplicados e fornecidos
pelos sistemas tratados.
24
2 Fundamentação Teórica
No próximo capítulo será apresentada e discutida a metodologia proposta e utilizada neste trabalho.
Capítulo 3
Seleção de estrutura através da análise de
resposta em frequência
3.1 Introdução
No capítulo anterior foram discutidas as etapas de Identif cação de Sistemas e a Transformada
de Fourier de sinais discretos.
Na seção 2.2.1.2 foi visto que as não linearidades de modelos NARX polinomiais são determinadas por produtos entre entradas, entre saídas ou por produtos entre entradas e saídas. Os
termos dos modelos são aqueles apresentados na equação (2.1). A escolha de estrutura destes
modelos consiste em determinar quais termos e agrupamentos farão parte do modelo, ou seja,
farão parte de sua estrutura com o objetivo de representar os regimes dinâmicos descritos pelos
dados.
É importante salientar que, se o número de termos do modelo for maior que o necessário, resultando em um modelo sobreparametrizado (modelo com um número excessivo de termos), seus
termos redundantes, além de acarretarem grande tempo computacional e mal condicionamento
numérico, podem fazer com que o modelo apresente regimes dinâmicos espúrios (AGUIRRE;
BILLINGS, 1995b; PIRODDI; SPINELLI, 2003) ou, até mesmo, torne-se instável. A escolha
de termos apropriados é de fundamental importância para obtenção de um modelo compacto e
robusto.
Alguns métodos tratados na literatura para seleção de estrutura foram apresentados na seção
2.2.1.3. Entre as técnicas usadas foram destacados o ERR que ordena os termos candidatos de
26
3 Seleção de estrutura através da análise de resposta em frequência
acordo com a contribuição de cada um na explicação da dinâmica dos dados. Tendo ordenado
os termos ainda é necessário achar o ponto de corte. Uma primeira estimativa pode ser obtida
através do critério de informação de Akaike. Como essas duas técnicas estão bem estabelecidas na literatura, serão chamadas de “métodos clássicos” daqui em diante para facilitar sua
referência no texto.
Neste capítulo será apresentada uma metodologia de suporte à seleção de estrutura de modelos
NARX polinomiais realizada a partir da análise no domínio da frequência de um sistema de uma
entrada e uma saída (SISO - do inglês Single Input, Single Output). O modelo deve reproduzir
aproximadamente a resposta em frequência do sistema em sua saída. Portanto, realizou-se um
estudo para saber como cada termo ou agrupamento deste colabora na reconstrução da resposta
em frequência do modelo que compõe. A análise matemática deste estudo é demonstrada na
seção 3.2.
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
O objetivo de se fazer a análise no domínio da frequência de um sistema é verif car como tal
sistema responde a um sinal senoidal de frequência única ω. Se for aplicado um sinal senoidal
de frequência única e média zero em um modelo polinomial NARX, esse modelo produzirá
frequências na saída que dependerão do conjunto de agrupamentos de termos que o compõem.
Nas próximas seções serão realizadas as análises matemáticas para os agrupamentos de termos
puros de entrada, agrupamentos de termos cruzados e agrupamento de termos puros de saída.
Ressalta-se que neste trabalho priorizou-se pela análise de ganhos de cada uma das frequências.
As análises de avanço e atraso de fase serão omitidas e sugeridas como assunto para trabalhos
futuros.
3.2.1 Agrupamento de termos puros de entrada
Considere um sinal de frequência única ω, média zero e ganho igual a um representado na forma
de Euler e apresentado na equação (3.1).
sen(kωT ) =
1 jkωT
(e
− e−jkωT ),
2j
(3.1)
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
27
sendo ω a frequência fundamental e T o período do sinal. Agora considere um modelo NARX
polinomial da equação (2.2) revisitado na equação (3.2).
y(k) =
y ,nu
m nX
X̀ X
m=0 p=0 n1 ,nm
cp,m−p (n1 , . . . ,nm )
p
Y
y(k − ni )
i=1
m
Y
u(k − ni ),
(3.2)
i=p+1
Aplicando o sinal senoidal apresentado na equação (3.1) em um modelo NARX como o apresentado na equação (3.2) pode-se observar o surgimento de frequências diferentes da frequência
de excitação do modelo. Tanto as frequências surgidas quanto o número dessas frequências dependerá de cada termo escolhido para a composição do modelo.
Q
Tratando-se de agrupamento de termos puros de entrada dados por c0,m (n1 , . . . ,nm ) m
i=1 u(k −
ni ) com p = 0, tem-se:
m
Y
1
c0,m (n1 , . . . ,nm ) ( (ej(k−ni )T − e−j(k−ni)T )).
2j
i=1
(3.3)
Sendo a equação (3.3) um termo puro de entrada de grau de não linearidade m, pode-se expandila utilizando o binômio de Newton, logo tem-se:
m
1 X
m!
·
(ejωT )m−g (−e−jωT )g .
2j g=0 g!(m − g)!
(3.4)
A equação (3.4) apresenta a expansão usando o binômio de Newton da equação (3.3) submetida
a um grau de não linearidade m. Essa equação mostra que dependendo de m um padrão de
frequências surge na saída. A partir dela é possível def nir exatamente quais frequências surgem
quando existe no modelo um agrupamento de termo puro de entrada. Os exemplos 3.2.1 a 3.2.4
apresentarão o desenvolvimento da expansão para agrupamento de termos puros de entrada de
graus de não linearidade ` variando de dois a cinco.
Para os exemplos a seguir considere o sinal de entrada apresentado na equação (3.1) na forma
de Euler.
Exemplo 3.2.1
Considere o termo θi u(k)2 de um modelo NARX qualquer.
Considerando o sinal aplicado como a entrada mostrada na equação (3.1) tem-se o desenvolvi-
28
3 Seleção de estrutura através da análise de resposta em frequência
mento:
1 jkωT
(e
− e−jkωT ))2
2j
1
= ( )2 · (ejkωT − e−jkωT )2
2j
1
= − · (e2jkωT + e−2jkωT − 2)
4
e2jkωT
e−2jkωT
1
= −
−
+
4
4
2
u(k)2 = (
(3.5)
Pela equação (3.5) pode-se notar que o agrupamento θi u(k)2 gera a frequencia zero e a frequência 2ω.
Exemplo 3.2.2
Agora considere o termo θi u(k)3 de um modelo NARX qualquer.
Considerando o sinal aplicado como a entrada mostrada na equação (3.1) tem-se o desenvolvimento:
1 jkωT
(e
− e−jkωT ))3
2j
1
= ( )3 · (ejkωT − e−jkωT )3
2j
1
= − · (−e−3jkωT − 3ejkωT + 3e−jkωT + e3jkωT )
8j
(e−3jkωT + 3ejkωT − 3e−jkωT − e3jkωT )
=
8j
u(k)3 = (
(3.6)
Pela equação (3.6) pode-se notar que o agrupamento θi u(k)3 gera a frequência ω e a frequência
3ω.
Exemplo 3.2.3
Agora considere o termo θi u(k)4 de um modelo NARX qualquer.
Considerando o sinal aplicado como a entrada mostrada na equação (3.1) tem-se o desenvolvi-
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
29
mento:
1 jkωT
(e
− e−jkωT ))4
2j
1
= ( )4 · (ejkωT − e−jkωT ))4
2j
1
· (6 − 4e2jkωT + e−4jkωT + e4jkωT − 4e−2jkωT )
=
16
(6 − 4e2jkωT + e−4jkωT + e4jkωT − 4e−2jkωT )
=
16
u(k)4 = (
(3.7)
Pela equação (3.6) pode-se notar que o agrupamento θi u(k)4 gera a frequência zero, a frequência
2ω e a frequência 4ω.
Exemplo 3.2.4
Considere o termo θi u(k)5 de um modelo NARX qualquer.
Considerando o sinal aplicado como a entrada mostrada na equação (3.1) tem-se o desenvolvimento:
1 jkωT
(e
− e−jkωT ))5
2j
1
= ( )5 · (ejkωT − e−jkωT )5
2j
1
=
· (−10e−jkωT + 10ejkωT − e−5jkωT − 5e3jkω + 5e−3jkωT + e5jkωT )
32j
(−10e−jkωT + 10ejkωT − e−5jkωT − 5e3jkωT + 5e−3jkωT T + e5jkωT )
=
(3.8)
32j
u(k)5 = (
Pela equação (3.6) pode-se notar que o agrupamento θi u(k)5 gera as frequências ω, 3ω e 5ω.
A partir dos desenvolvimentos dos exemplos 3.2.1 a 3.2.4 pôde ser visto que agrupamentos
de termos puros de entrada de grau de não linearidade ` par não geram ω. Apenas geram
frequências múltiplas pares de ω. Os mesmos agrupamentos, porém de grau de não linearidade
` ímpar geram a frequência ω e mais frequências múltiplas ímpares de ω.
Um exemplo gráf co será apresentado na Figura 3.1. Neste exemplo são apresentados as análises
30
3 Seleção de estrutura através da análise de resposta em frequência
no domínio da frequência fornecidas por cada agrupamento de termos puros de entrada demonstrados nos exemplos 3.2.1 a 3.2.4.
Exemplo 3.2.5
Considere os agrupamentos de termos puros de entrada apresentados nos exemplos 3.2.1 a 3.2.4.
Aplicando-se um sinal de entrada senoidal de média zero amostrado a uma taxa de 62,8319
rad/s de ganho um e frequência única de 1 rad/s, obteve-se a DTFT para cada um deles. Os
coef cientes de todos os termos representados foi def nido como um. As análises no domínio da
frequência para os exemplos são apresentadas na Figura 3.1.
Na Figura 3.1 (b), é apresentada a resposta em frequência de um termo linear de entrada. Podese observar que a frequência fornecida pelo termo é a mesma frequência aplicada na entrada
apresentada na Figura 3.1 (a). O mesmo não se pode dizer das respostas em frequência exibidas nas Figuras 3.1 (c) a (f). A Figura 3.1 (c) refere-se a uma resposta em frequência de
um agrupamento de termos puro de entrada não-linear de grau de não linearidade ` igual a dois.
Considerando que agrupamentos com ` par produzem as frequências kω sendo k = 2,4,6, . . . ,m
na saída e a frequência fundamental é de 1 rad/s, logo este agrupamento fornecerá a frequência
2 rad/s.
A mesma regra é válida para o caso da Figura 3.1 (d) onde é apresentada a resposta em frequência de um agrupamento em que ` é igual a quatro. Neste caso surgem as duas frequências:
zero, 2 rad/s e 4 rad/s. Tanto para o caso de ` é igual a dois e ` é igual a quatro, as saídas não
reproduzem a frequência fundamental.
Para os casos representados pelas Figuras 3.1 (e) e (f) os agrupamentos de termos puros de entrada não lineares apresentam grau de não linearidade ` ímpar. Na Figura 3.1 (e) é apresentada a
resposta em frequência de um agrupamento de ` igual a três. Na Figura 3.1 (f), ` é igual a cinco.
Como agrupamentos puros de entrada de não linearidade ` ímpar produzem as frequências kω
sendo k = 1,3,5, . . . ,m na saída, serão produzidas pelo agrupamento representado na Figura
3.1 (e) a frequência fundamental ω e a frequência 3ω, e para o agrupamento representado na
Figura 3.1 (f), a frequência fundamental ω e as frequências 3ω e 5ω.
A partir da discussão e dos exemplos 3.2.1 a 3.2.5 pode-se generalizar: para agrupamentos de
termos puros de entrada de grau de não linearidade m:
• sendo m par, serão produzidas as frequências 2kω, sendo k = 0,1 . . . ,m/2, m ≥ 2 e;
• sendo m ímpar, serão produzidas as frequências 2(k + 1)ω, sendo k = 0,1 . . . ,(m − 1)/2,
m ≥ 3.
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
(a)
(b)
Espectro de frequência da saída (y)
0.8
0.8
0.7
0.7
0.6
0.6
0.5
0.5
|Y(jω)|
|U(jω)|
Espectro de frequência da entrada (u)
0.4
0.4
0.3
0.3
0.2
0.2
0.1
0.1
0
0
1
2
3
Frequencia (rad/s)
4
5
0
6
0
1
2
(c)
4
5
6
5
6
5
6
Espectro de frequência da saída (y)
0.45
0.45
0.4
0.4
0.35
0.35
0.3
|Y(jω)|
0.3
|Y(jω)|
3
Frequencia (rad/s)
(d)
Espectro de frequência da saída (y)
0.25
0.25
0.2
0.2
0.15
0.15
0.1
0.1
0.05
0.05
0
31
0
1
2
3
Frequencia (rad/s)
4
5
0
6
0
1
2
(e)
3
Frequencia (rad/s)
4
(f)
Espectro de frequência da saída (y)
Espectro de frequência da saída (y)
0.6
0.5
0.45
0.5
0.4
0.35
|Y(jω)|
|Y(jω)|
0.4
0.3
0.3
0.25
0.2
0.2
0.15
0.1
0.1
0.05
0
0
1
2
3
Frequencia (rad/s)
4
5
6
0
0
1
2
3
Frequencia (rad/s)
4
Figura 3.1: Análise no domínio da frequência de agrupamentos de termos de entrada submetidos a um
sinal de frequência única. Em (a) a DTFT da entrada é apresentada e em (b) a DTFT do
sinal fornecido por um termo linear de entrada. Em (c) e em (d) são apresentadas as análises
no domínio da frequência de sinais gerados por agrupamento de termos com graus de não
linearidade pares sendo ` = 2 para (c) e ` = 4 para (d). Finalmente, em (e) e em (f) são
apresentadas as análises no domínio da frequência de sinais gerados por agrupamento de
termos com graus de não linearidade ímpares sendo ` = 3 para (e) e ` = 5 para (f))
32
3 Seleção de estrutura através da análise de resposta em frequência
Portanto, foram def nidas quais as frequências que um agrupamento de termos puros de entrada
fornecem para a resposta do modelo quando submetidos a uma sinal de frequência única. Graus
de não linearidade ` > 5 podem ser analisados expandindo-se a equação (3.4) comprovando a
generalização proposta.
3.2.2 Agrupamento de termos cruzados
Considerando o modelo NARX polinomial da equação (3.2), pode-se notar que agrupamento de
termos cruzados de entrada e saída também podem compor o modelo. Representado este tipo
de agrupamento como na equação (3.9), tem-se:
cp,m−p (n1 , . . . ,nm )
p
Y
i=1
y(k − ni )
m
Y
u(k − ni )
(3.9)
i=p+1
com p,m 6= 0.
Considere o exemplo 3.2.6.
Exemplo 3.2.6
Sendo a equação (3.10) parte de um modelo NARX polinomial de grau de não linearidade `
igual a dois estimado e com parâmetros θ̂ tem-se:
y(k) = θ̂i u(k)y(k − 1) + θ̂i+1 u(k).
(3.10)
Aplicando-se o mesmo sinal de entrada apresentado na equação (3.10) e considerando que y(k−
1) = 0, tem-se:
y(k) = θ̂i (
no passo k, e:
1 jkωT
1
(e
− e−jkωT )) · 0 + θ̂i+1 ( (ejkωT − e−jkωT )))
2j
2j
(3.11)
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
33
1 j(k+1)ωT
1
(e
− e−j(k+1)ωT )) · (θi ( (ej(k+1)ωT − e−j(k+1)ωT )) · 0 +
2j
2j
1 j(k+1)ωT
+θi+1 ( (e
− e−j(k+1)ωT )) +
2j
1
+θi+1 ( (ej(k+1)ωT − e−j(k+1)ωT ))
2j
1
1
= θi ( (ej(k+1)ωT − e−j(k+1)ωT )) · (θi+1 ( (ej(k+1)ωT − e−j(k+1)ωT ))) +
2j
2j
1 j(k+1)ωT
+θi+1 ( (e
− e−j(k+1)ωT ))
2j
1
1
= θi θi+1 ( (ej(k+1)ωT − e−j(k+1)ωT ))( (ej(k+1)ωT − e−j(k+1)ωT )) +
2j
2j
1 j(k+1)ωT
− e−j(k+1)ωT ))
+θi+1 ( (e
2j
e2j(k+1)ωT
e−2j(k+1)ωT
1
θi θi+1 (−
−
− )+
4
4
2
1 j(k+1)ωT
−j(k+1)ωT
= +θi+1 ( (e
−e
))
(3.12)
2j
y(k + 1) = θi (
no passo k + 1, e
e2j(k+2)ωT
e−2j(k+2)ωT
1
1 j(k+2)ωT
(e
− e−j(k+2)ωT )) · (θi θi+1 (−
−
− )+
2j
4
4
2
1
1
+θi+1 ( (ej(k+2)ωT − e−j(k+2)ωT )) + θi+1 ( (ej(k+2)ωT − e−j(k+2)ωT ))
2j
2j
−3j(k+2)ωT
j(k+2)ωT
−j(k+2)ωT
(e
+ 3e
− 3e
− e3j(k+2)ωT )
= θi2 θi+1 (
)
8j
e−2j(k+2)ωT
1
e2j(k+2)ωT
−
− )+
+θi θi+1 (−
4
4
2
1 j(k+2)ωT
+θi+1 ( (e
− e−j(k+2)ωT ))
(3.13)
2j
y(k + 2) = θi (
no passo k + 2.
Pode-se notar através dos desenvolvimentos apresentados nas equações (3.11) a (3.13), que a
cada iteração o modelo é realimentado por uma saída passada que já contém frequências in-
34
3 Seleção de estrutura através da análise de resposta em frequência
seridas pelo próprio modelo em passos anteriores. Isso faz com que o sinal gerado por um
agrupamento de termos cruzados forneça à própria saída inf nitas frequências ressonantes da
frequência fundamental. Quanto maior o número de passos dados pelo modelo, mais frequências surgirão na análise no domínio da frequência.
Exemplo 3.2.7
Neste exemplo é apresentado a análise no domínio da frequência de um agrupamento de termos
cruzados. Aplicando-se um sinal de entrada senoidal de média zero amostrado a uma taxa
de 62,8319 rad/s de ganho um e frequência única de 1 rad/s, obteve-se a DTFT da saída. O
coef ciente do agrupamento do exemplo é um. As análises no domínio da frequência para o
exemplo 3.2.6 são apresentadas na Figura 3.2.
(a)
(b)
Espectro de frequência da entrada (u)
Espectro de frequência da saída (y)
3
0.8
0.7
2.5
0.6
|Y(jω)|
|U(jω)|
2
0.5
0.4
0.3
1.5
1
0.2
0.5
0.1
0
0
1
2
3
4
5
Frequencia (rad/s)
6
7
8
0
0
1
2
3
4
5
Frequencia (rad/s)
6
7
8
Figura 3.2: Análise no domínio da em frequência de um modelo com agrupamentos de termos cruzados
submetidos a um sinal de frequência única. Em (a) a DTFT da entrada é apresentada e em
(b) a DTFT do sinal fornecido por um agrupamento de termos cruzados.
Pela Figura 3.2 é possível notar que aplicando uma frequência na entrada, várias outras frequências surgem na saída. Isso pode ser comprovado se o desenvolvimento das equações (3.11)
a (3.13) for continuado.
Portanto, para agrupamento de termos cruzados de grau de não linearidade m, pode-se generalizar, a partir dos exemplos 3.2.6 e 3.2.7, que:
• qualquer que seja m, o número de frequências fornecido por um agrupamento de termos
cruzados, aumenta à medida que as iterações com saídas passadas cresce. Logo, quando
k → inf, são fornecidas inf nitas frequências na saída múltiplas de ω.
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
35
3.2.3 Agrupamento de termos puros de saída
Para agrupamento de termos puros de saída a reprodução de frequências na saída é semelhante
à reprodução de frequências para os agrupamentos de termos cruzados.
Considerando novamente o modelo NARX polinomial da equação (3.2). Representando o agrupamento de termos puros de saída como na equação (3.9), tem-se:
cp,0(n1 , . . . ,np )
p
Y
y(k − ni )
(3.14)
i=1
com p 6= 0.
O exemplo 3.2.8 apresentará o desenvolvimento de um modelo qualquer que contém um agrupamento de termos puro de saída. O grau de não linearidade do termo do exemplo é dois mas
pode ser expandido para graus de não linearidade pares superiores. Para o exemplo 3.2.8 serão
considerados graus de não linearidades pares para os agrupamentos.
Exemplo 3.2.8
Sendo a equação (3.15) parte de um modelo NARX de grau de não linearidade par ` igual a dois
estimado e com parâmetros θ̂ tem-se:
ŷ(k) = θ̂i y(k − 1)2 + θ̂i+1 u(k).
(3.15)
Aplicando-se o mesmo sinal de entrada apresentado na equação (3.15) e considerando que y(k−
1) = 0, tem-se:
ŷ(k) = θ̂i · (0)2 + θ̂i+1 (
no passo k, e:
1 jkωT
(e
− e−jkωT )))
2j
(3.16)
36
3 Seleção de estrutura através da análise de resposta em frequência
1 j(k+1)ωT
1
(e
− e−j(k+1)ωT )))2 + θ̂i+1 ( (ej(k+1)ωT − e−j(k+1)ωT ))
2j
2j
−2j(k+1)ωT
2j(k+1)ωT
e
1
1
e
2
−
− ) + θ̂i+1 ( (ej(k+1)ωT − e−j(k+1)ωT ))
= θ̂i θ̂i+1
(−
4
4
2
2j
ŷ(k + 1) = θ̂i (θ̂i+1 (
no passo k + 1, e:
2
ŷ(k + 2) = θ̂i (θ̂i θ̂i+1
(−
e2j(k+2)ωT
e−2j(k+2)ωT
1
1
−
− ) + θ̂i+1 ( (ej(k+2)ωT − e−j(k+2)ω )))2
4
4
2
2j
1 j(k+2)ωT
(e
− e−j(k+2)ωT ))
2j
e2j(k+2)ωT
e−4j(k+2)ωT
e4j(k+2)ωT
3
e−2j(k+2)ωT
3 4
+
+
+
+ )) +
= θ̂i θ̂i+1 (
4
4
16
16
8
−j(k+2)ωT
j(k+2)ωT
−3j(k+2)ωT
j(k+2)ωT
e
e
j
e
j
e
j
3
+2θ̂i2 θ̂i+1
(
−
+
+
−
4j
4
8
8
e−j(k+2)ωT
e3j(k+2)ωT
1 e2j(k+2)ωT
e−2j(k+2)ωT
2
−
−
) + θ̂i θ̂i+1
( −
−
)+
8j
8
2
4
4
1
+θ̂i+1 ( (ej(k+2)ωT − e−j(k+2)ωT ))
2j
+θ̂i+1 (
(3.17)
no passo k + 2.
Nota-se que há o surgimento de todas as frequências ressonantes da frequência fundamental.
Se houver continuidade no desenvolvimento para os passos posteriores, haverá o surgimento
de inf nitas frequências quando um agrupamento de termos puro de saída com grau de não
linearidade dois compor o modelo.
No exemplo 3.2.9 será demonstrado o agrupamento de termos puro de saída para grau de não
linearidade ` = 2 ímpar igual a três. Tal processo pode ser estendido para qualquer grau de não
linearidade ímpar para agrupamentos puros de saída.
Exemplo 3.2.9
Sendo a equação (3.18) parte de um modelo NARX de grau de não linearidade ` igual a dois
estimado e com parâmetros θ̂ tem-se:
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
ŷ(k) = θ̂i y(k − 1)3 + θ̂i+1 u(k).
37
(3.18)
Aplicando-se o mesmo sinal de entrada apresentado na equação (3.15) e considerando que y(k−
1) = 0, tem-se:
ŷ(k) = θ̂i · (0)3 + θ̂i+1 (
1 jkωT
(e
− e−jkωT )))
2j
(3.19)
no passo k, e:
1 j(k+1)ωT
1
(e
− e−j(k+1)ωT )))3 ) + θ̂i+1 ( (ej(k+1)ωT − e−j(k+1)ωT )))
2j
2j
1
3
= θ̂i θ̂i+1
(− · (−e−3j(k+1)ωT − 3ej(k+1)ωT + 3e−j(k+1)ωT + e3j(k+1)ωT ))
8j
1 j(k+1)ωT
+θ̂i+1 ( (e
− e−j(k+1)ωT )))
(3.20)
2j
ŷ(k + 1) = θ̂i (θ̂i+1 (
no passo k + 1, e:
3
ŷ(k + 2) = θ̂i (θ̂i θ̂i+1
(−
1
· (−e−3j(k+2)ωT − 3ej(k+2)ωT + 3e−j(k+2)ωT + e3j(k+2)ωT ))
8j
1 j(k+2)ωT
(e
− e−j(k+2)ωT ))))3 +
2j
1
+θ̂i+1 ( (ej(k+2)ωT − e−j(k+2)ωT ))
2j
+θ̂i+1 (
(3.21)
no passo k + 2. A expansão da parte da equação (3.21) elevada ao cubo gera uma outra equação
com muitos termos que pode ser expandida computacionalmente e será suprimida aqui. O que é
relevante é que essa expansão gerará frequências múltiplas ímpares da frequência fundamental.
Exemplo 3.2.10
Neste exemplo é apresentado a análise no domínio da frequência de um agrupamento de termos
38
3 Seleção de estrutura através da análise de resposta em frequência
puros de saída. Aplicando-se um sinal de entrada senoidal de média zero amostrado a uma taxa
de 62,8319 rad/s de ganho um e frequência única de 1 rad/s, obteve-se a DTFT da saída. O
coef ciente do agrupamento do exemplo é 0,1. As análises no domínio da frequência para o
exemplo 3.2.6 são apresentadas na Figura 3.3.
(a)
(b)
Espectro de frequência da entrada (u)
Espectro de frequência da saída (y)
0.09
0.7
0.08
0.6
0.07
0.06
0.5
|Y(jω)|
|U(jω)|
0.1
0.8
0.4
0.05
0.04
0.3
0.03
0.2
0.02
0.1
0
0.01
0
1
2
3
4
5
Frequencia (rad/s)
6
7
0
8
0
1
2
3
(c)
0.1
0.09
0.04
0.08
0.035
0.07
0.03
0.06
|Y(jω)|
|Y(jω)|
0.05
0.025
0.04
0.015
0.03
0.01
0.02
0.005
0.01
2
3
4
5
Frequencia (rad/s)
8
9
7
8
9
0.05
0.02
1
7
Espectro de frequência da saída (y)
0.045
0
6
(d)
Espectro de frequência da saída (y)
0
4
5
Frequencia (rad/s)
6
7
8
9
0
0
1
2
3
4
5
Frequencia (rad/s)
6
Figura 3.3: Análise no domínio da frequência de um modelo com agrupamentos de termos cruzados
submetidos a um sinal de frequência única. Em (a) a DTFT da entrada é apresentada. Em
(b) a DTFT do sinal fornecido por um agrupamento de termos puro de saída com grau de
não linearidade quatro. Em (c) a DTFT do sinal fornecido por um agrupamento de termos
puro de saída com grau de não linearidade três. E em (d) a DTFT do sinal fornecido por um
agrupamento de termos puro de saída com grau de não linearidade cinco.
Pela Figura 3.3 (b) é possível notar que aplicando uma frequência na entrada, várias outras
frequências surgem na saída para um grau de não linearidade par igual a quatro em um agrupamento de termos puro de saída. Desenvolvimentos e simulações com outros graus de não
3.2 Análise matemática da resposta em frequência de modelos NARX polinomiais
39
linearidades pares comprovam este comportamento.
Na Figura 3.3 (c) e (d) pode-se notar que apenas frequências múltiplas ímpares da frequência
fundamental são geradas além da frequência fundamental. Tanto para o grau de não linearidade
do agrupamento igual a três como para o grau de não linearidade igual a cinco. Nestes casos o
que diferencia o caso de grau de não linearidade três para o cinco é o ganho das frequências. O
aumento do grau de não linearidade ref ete num aumento do ganho das frequências ressonantes
da frequências fundamental reproduzidas pelo agrupamento.
Portanto, para agrupamento de termos puros de saída de grau de não linearidade m, pode-se
generalizar, a partir dos exemplos 3.2.8, 3.2.9 e 3.2.10, que:
• para graus de não linearidade ` pares são geradas inf nitas frequências múltiplas da frequência fundamental e;
• para graus de não linearidade ` ímpar são geradas inf nitas frequências múltiplas ímpares
da frequência fundamental.
3.2.4 Considerações sobre a análise matemática no domínio da frequência de modelos
NARX polinomiais
Nos exemplos apresentados na seção anterior alguns modelos foram simulados para verif car a
análise da resposta em frequência para cada tipo de agrupamento de termos individualmente.
Para que o modelo convergisse, apenas termos lineares de entrada foram adicionados além do
agrupamento em estudo. Logo, não houve interferência nas frequências produzidas para cada
agrupamento.
Pôde ser visto que agrupamentos de termos não lineares puros de entrada submetidos a um
sinal de entrada de frequência única têm em sua saída frequências bem def nidas. Pela análise
matemática é possível notar que é simples também obter os ganhos e as defasagens.
No entanto, não se pode af rmar o mesmo para agrupamentos de termos cruzados. Pôde-se
comprovar que a cada iteração esses agrupamentos fornecem uma nova frequência ao sinal de
saída do modelo. Portanto, f ca difícil def nir quantas e quais frequências serão fornecidas por
cada um destes agrupamentos.
Para agrupamentos de termos puros de saída foi visto que o número de frequências geradas na
saída não pode ser def nido. No entanto, pôde-se observar que esses agrupamentos com grau de
não linearidade ímpar fornecem a frequência fundamental a suas frequências múltiplas ímpares
ressonantes.
40
3 Seleção de estrutura através da análise de resposta em frequência
É importante notar também que os ganhos de cada frequência dependem dos coef cientes que
acompanham cada agrupamento. Agrupamentos diferentes podem fazer surgir as mesmas frequências na saída do modelo, logo os ganhos das frequências apresentadas é resultado de uma
soma ponderada pelos coef cientes de cada termo que gera a mesma frequência.
Uma outra observação que deve ser feita é que os atrasos de cada termo do agrupamento não
interferem em quais frequências surgem na saída devido à presença do agrupamento. Isso pode
ser comprovado pelos desenvolvimentos realizados nesta seção. A Tabela 3.1 apresenta a relação de termos e o surgimento de frequências.
Tabela 3.1: Resposta em frequência de agrupamentos de termos submetidos a um sinal de entrada de
frequência única ωi .
Agrupamento
Frequências de saída
Ωu
ωi
Ωy
ωi
Ωum m par
kωi sendo que k = 0,2,4,6, . . . ,m
Ωyp p par
kωi sendo que k = 1,2,3, . . . ,∞
Ωum m ímpar kωi sendo que k = 1,3,5, . . . ,m
Ωyp p ímpar
kωi sendo que k = 1,3,5, . . . ,∞
Ωum yp
kωi sendo que k = 1,2,3, . . . ,∞
Conforme apresentado na Tabela 3.1 pode-se notar que def nir as frequências surgidas na saída
para cada agrupamento não é trivial. Apenas para os termos lineares e para os agrupamentos de
termos puros de entrada é possível chegar a um número def nido de frequências geradas.
Quando trata-se de ganho de cada uma dessas frequências surgidas, é necessário conhecer o
coef ciente que acompanha cada uma das fontes dessas frequências, no caso, cada um dos agrupamentos.
Se um mesmo modelo contiver termos lineares, agrupamentos de termos puros de entrada com
grau de não linearidade ímpar, agrupamento de termos cruzados e agrupamento de termos puros
de saída, todos esses irão contribuir para o ganho da frequência fundamental ωi . Quanto às
frequências ressonantes, os agrupamentos de termos puros de entrada irão contribuir para as
frequências ímpares múltiplas de ωi e os demais agrupamentos para todas as outras. Se contiver
agrupamentos de termos puros de entrada com grau de não linearidade par, este contribuirá para
os ganhos das frequências múltiplas de 2ωi.
Como um modelo deve aproximar a resposta em frequência do sistema modelado, a estimação
de parâmetros deve ser realizada de forma a ponderar os ganhos das frequências geradas por
3.3 Metodologia aplicada
41
cada agrupamento e esta soma f nal aproximar ao máximo da frequência real de saída gerada
pelo sistema. O mesmo deve ser considerado para os avanços e/ou atrasos de fase.
Como se pode notar pelos desenvolvimentos realizados ao longo da seção 3.2.1, um modelo já
parametrizado, será capaz de fornecer apenas um padrão de resposta em frequência para cada
frequência individualmente aplicada na entrada. Logo, se for aplicado um sinal de frequência
única na entrada de um modelo, este sempre fornecerá na saída o mesmo padrão de reposta em
frequência. A partir dessa informação é proposta a metodologia utilizada neste trabalho.
3.3 Metodologia aplicada
A metodologia aplicada neste trabalho baseia-se nos desenvolvimentos realizados na seção 3.2
e serão considerados apenas o ganho das frequências surgidas nas saídas desconsiderando-se
o avanço e/ou atraso de fase que possam ocorrer. Pode-se observar que um modelo NARX
polinomial tem o mesmo padrão de resposta em frequência quando é submetido a uma única
frequência. Isso signif ca que o número de frequências ressonantes que surgirão na saída do
modelo será sempre o mesmo para qualquer frequência que for aplicada individualmente. Isso
também vale para o ganho de cada uma dessas frequências desde que não se localizem em uma
região de atenuação.
Entretanto, um sistema pode f ltrar determinadas frequências. Se sua dinâmica for lenta, pode
não haver resposta a sinais com componentes de frequências muito altas. Portanto, certas
regiões de frequências de excitação podem não interferir na resposta do sistema. Algumas frequências podem ser atenuadas pelo próprio sistema e cabe ao modelo reproduzir esta atenuação
em sua saída.
Logo, um modelo NARX que represente bem um sistema, além de gerar as frequências ressonantes na saída, precisa delimitar a sua região de frequências de interesse. Essa delimitação é
realizada pelo conjunto de termos lineares que formam um tipo de filtro que funciona como um
bloco funcional para suprimir frequências espúrias dos sistemas nos modelos.
3.3.1 Sugestão de grau de não linearidade ` a partir da análise no domínio da frequência
Um sistema que reproduza na saída, apenas a mesma frequência aplicada na entrada pode ser
representado por um modelo de grau de não linearidade ` igual a um, ou seja, linear. Para
sistemas que reproduzam na saída frequências ressonantes da frequência fundamental, pode-se
fazer uma análise partindo da discussão realizada na seção 3.2.1 para agrupamentos de termos
42
3 Seleção de estrutura através da análise de resposta em frequência
puros de entrada.
Se um sistema apresenta a frequência fundamental mais uma frequência ressonante pode-se
sugerir, por exemplo, que o grau de não linearidade seja dois e que se tenha termos de entrada
lineares que serão responsáveis pela reprodução da frequência fundamental. Caso o sistema
apresente duas frequências ressonantes pode-se sugerir um grau de não linearidade três e pelo
menos um termo com grau de não linearidade dois. O termo com grau de não linearidade três
faz surgir na saída a frequência fundamental e sua primeira frequência múltipla ímpar. O termo
com grau de não linearidade par faz surgir a frequência par múltipla da frequência fundamental.
Ainda no caso de haver duas frequências ressonantes, poderiam ser sugeridos agrupamentos de
termos puros de entrada de grau de não linearidade dois e agrupamentos de termos cruzados de
entrada e saída visto que estes são capazes de gerar todas as frequências ressonantes. Adicionar
agrupamentos de termos cruzados ao modelo pode também colaborar na representação de frequências ressonantes, podendo-se assim, reduzir a complexidade do modelo no que diz respeito
ao grau de não linearidade. No entanto, pode-se fazer surgir frequências espúrias ao sistema
representado.
O grau de não linearidade do modelo pode partir da análise no domínio da frequência dos agrupamentos de termos de entrada. Isso garantirá que todas as frequências dentro da região de
potência espectral sejam representadas. Como foi visto na seção 3.2.1, os agrupamentos de
termos puros de entrada geram frequências que podem ser def nidas. Tomando as informações
da Tabela 3.1 e comparando com a análise no domínio da frequência do sistema, é possível
determinar o grau de não linearidade do agrupamento de termos puros de entradas que represente todas as frequências apresentadas. Então toma-se este grau de não linearidade como
máximo garantindo-se a representação de todas as frequências da saída e limita-se o grau de
não linearidade dos demais agrupamentos de termos a este máximo.
De acordo com a Tabela 3.1 pode-se notar que a simples presença de agrupamentos de termos
cruzados contribui para o aparecimento de diversas frequências na saída quando o modelo é
submetido a um sinal de entrada de frequência única. Logo, pode-se pressupor que um grau de
não linearidade dois ou três em um modelo que contenha agrupamento de termos cruzados já é
possível representar todas as frequências ressonantes. Embora isso possa acarretar o surgimento
de frequências espúrias ao sistema modelado. Diante disso, uma boa representação já pode ser
obtida a partir de um grau de não linearidade dois ou três mesmo apresentando mais que duas
frequências ressonantes.
Porém, considerar o grau de não linearidade do sistema a partir dos agrupamentos de termos
puros de entrada é mais conveniente, pois, pode-se determinar exatamente quantas e quais fre-
3.3 Metodologia aplicada
43
quências são inseridas na saída do modelo quando é inserido um agrupamento deste tipo.
De uma forma generalizada, pode-se adotar que a quantidade de frequências fundamental e
ressonantes de um sistema surgidas em sua região de potência espectral pode determinar o grau
de não linearidade de um modelo NARX polinomial que representará bem o sistema.
Uma vez que ωs = 2ωM , logo, ωM = 1/2ωs , ou seja, 50%. No entanto, as taxas de amostragem
adotas foram levemente aumentadas neste trabalho. Portanto, considerou-se neste trabalho as
frequências ressonantes até o valor de 30% da frequência de amostragem. Logo, ` = nf r + 1
sendo nf r o número de frequências ressonantes até 30% da frequência de amostragem do sinal.
3.3.2 Seleção de estrutura a partir da resposta em frequência da saída do sistema
Um modelo que represente bem o sistema precisa reproduzir adequadamente a sua resposta em
frequência. No entanto, avaliar um sistema submetido a um sinal com diversas frequências na
entrada é inviável matematicamente. Além disso, pode-se notar através dos desenvolvimentos
e discussões realizados nas seções 3.2.1 a 3.2.3 que um modelo polinomial NARX apresentará
sempre os mesmos padrões de ganhos e provavelmente de avanço/atraso de fase para qualquer
que seja a frequência aplicada individualmente em sua entrada. Considerando que um sinal
de múltiplas frequências seja a soma de vários sinais de frequência única, logo a saída de um
modelo submetido a um sinal de múltiplas frequências será a soma no domínio da frequência
das saídas obtidas correspondentes a cada sinal de frequência única.
Portando, optou-se por realizar o procedimento proposto neste trabalho aplicando-se um sinal de
entrada de frequência única, observar a análise no domínio da frequência do sistema e comparála ao domínio da frequência do modelo ao acrescentar cada termo. Logo, um dos pré-requisitos
para o uso do método proposto é que o sistema a ser modelado permita aplicar um sinal de
frequência única em sua entrada com uma frequência próxima a menor frequência de potência
espectral do sistema e que excursione toda a faixa de entrada em amplitude.
Resumidamente, a metodologia proposta é a de realizar uma comparação entre as frequências geradas pelo sistema em sua saída e escolher os agrupamentos de termos que devem estar
presentes no modelo a partir da resposta a um sinal de frequência única. O método oferece
termos ao modelo e verif ca se houve uma aproximação do domínio da frequência do modelo
em relação ao do sistema. Neste trabalho, quando se refere ao modelo aproximar o domínio da
frequência do sistema, limita-se a análise do ganho de cada frequência.
A metodologia adotada para o desenvolvimento desse procedimento é sintetizada por etapas,
apresentadas a seguir e logo em seguida serão discutidas cada uma das etapas:
44
3 Seleção de estrutura através da análise de resposta em frequência
1. Realize as coletas de dados do sistema para identif cação;
(a) Separe dados para identif cação e para validação;
2. Obtenha a razão entre análise no domínio da frequência para delimitar a região de potência espectral do sistema;
3. Submeta o sistema a um sinal de frequência única entre 10 a 15 vezes menor que a maior
frequência de potência espectral;
4. Obtenha a análise no domínio da frequência do sistema;
5. PARA - Determine o grau de não linearidade (`) conforme o número de frequências
ressonantes produzidas na análise no domínio da frequência do passo anterior até o limite
de 30% da frequência de amostragem;
(a) Determine os máximos atrasos de entrada (nu ) e de saída (ny ) utilizando algum
método presente na literatura;
(b) Hierarquize os agrupamentos de termos utilizando o ERR;
(c) Sugira o acréscimo ou retirada dos agrupamentos de termos conforme discutido na
seção 3.2.4;
(d) Inicialize este laço com os termos lineares;
(e) ENQUANTO -
ErroW (i+1)−ErroW (i)
ErroW (i)
× 100 < 1% E ErroW (i + 1) < ErroW (i)
i. Adicione o próximo agrupamento de termos mais importante segundo o ERR;
ii. Realize a estimação de parâmetros;
iii. Simule o modelo submetido ao mesmo sinal de frequência única escolhido no
passo 3;
iv. Obtenha a análise no domínio da frequência do modelo submetido ao mesmo
sinal de frequência única escolhido no passo 3;
v. Separe as frequências fundamental e ressonantes até o limite de 30% da frequência de amostragem;
vi. Calcule ErroW usando a equação (3.22) para as frequências escolhidas no
passo anterior;

2
Grfi −Gmfi n+1
X Grfi ErroW = 
× Grfi 
Pn+1
Grf
j
j=1
i=1
(3.22)
3.3 Metodologia aplicada
45
sendo Grf o ganho de cada frequência do sistema a ser modelado, Gmf o
ganho de cada frequência do modelo da saída submetidos ao sinal de entrada
def nido no passo 3 e n o número de frequências consideradas.
FIM ENQUANTO
(f) Simule o modelo de estrutura def nida no passo anterior;
6. Caso o desempenho do modelo não atenda;
(a) Adicione 1 ao grau de não linearidade;
(b) Volte ao passo 5a;
7. Valide o modelo obtido (Seção 2.2.1.7).
Apresentado o algoritmo que sintetiza o procedimento proposto neste trabalho, são necessários
alguns comentários.
Ao realizar a coleta de dados, é necessário seguir o procedimento descrito neste trabalho na
seção 2.2.1.1. É importante observar que o sinal aplicado na entrada tem que conter o máximo
de frequências possível e se aproximar de uma entrada branca. Isso permitirá que se delimite
com maior exatidão as frequências de interesse do sistema.
Para obter as razões entre os ganhos das frequências vistas na análise no domínio da frequência
optou-se por utilizar a medida expressa na equação:
GdB = 20log10
|fy (jω)|
|fu (jω)|
(3.23)
sendo fy (jω) o vetor de frequências da saída gerada pelo modelo e fu (jω) é o vetor de frequências da entrada aplicada ao modelo obtidos através da DTFT.
Seguindo para o próximo passo, para determinar os máximos atrasos dos sinais de entrada
e saída do modelo, é necessário seguir algum procedimento presente na literatura. Um bom
exemplo é o método dos falsos vizinhos (FNN) onde os máximos atrasos são estimados diretamente a partir dos dados (RHODES; MORARI, 1995a,b) ou as funções de correlação cruzada
(FCC) de ordem elevada presentes em (BILLINGS; ZHU, 1994; AGUIRRE, 1997) e utilizados
neste trabalho.
A sugestão de termos do passo 5c pode ser feita seguindo as análises feitas na seção 3.2.4. A
partir do número de frequências ressonantes é possível determinar o grau de não linearidade de
partida para a representação segundo a discussão da seção 3.3.1.
46
3 Seleção de estrutura através da análise de resposta em frequência
As frequências que devem ser consideradas na saída tanto para a sugestão do grau de não linearidade quanto para comparação realizada no método são aquelas que estão na região de potência
espectral do sistema. Porém, é necessário avaliar quais frequências realmente fazem parte da
dinâmica do sistema e quais frequências surgem em função da amostragem e do janelamento
do sinal. Frequências que surgem após este limite provavelmente são espúrias à dinâmica do
sistema. Além disso deve-se optar por uma frequência o mais próximo da menor frequência
de potência espectral possível. Isso permitirá a análise de todas as frequências ressonantes à
frequência aplicada pelo sinal de entrada.
Pelos desenvolvimentos realizados nas seções 3.2.1 a 3.2.3 e resumidos na Tabela 3.1 podem
ser sugeridas com maior exatidão apenas os agrupamentos de termos puros de entrada.
Pode-se também, ao invés de adotar o procedimento proposto como ferramenta principal na
seleção de estrutura, utilizá-lo como ferramenta auxiliar a outros métodos onde estruturas com
menores números de termos já foram escolhidas cabendo ao método realizar o corte de termos
menos relevantes.
Uma boa prática nesse caso, seria utilizar o método como auxílio aos critérios ERR e Akaike.
Utilizando esses critérios antes da aplicação da ferramenta proposta o método avalia se há
relevância de todos os termos incluídos nas estruturas pelos métodos clássicos. Portanto, seria uma opção ao método acrescentar o seguinte passo antes do passo 5c:
• Utilize o critério de corte de Akaike.
Portanto, o método iria buscar por estruturas mais compactas do que aquelas oferecidas pela
aplicação do critério ERR em conjunto com o critério de corte de Akaike.
A estimação de parâmetros do passo 5(e)ii pode ser feita utilizando o MQ ou MQE discutidos
na seção 2.2.1.6.
É importante apontar agrupamentos de termos importantes na explicação da dinâmica do sistema modelado porque diminuirá a variedade de estruturas avaliadas no passo 5e diminuindo o
tempo de processamento computacional.
O modelo agora com os parâmetros estimados precisa reproduzir o mesmo domínio da frequência que o sistema submetido ao mesmo sinal de entrada. Para isso, a análise no domínio da
frequência é obtido utilizando a DTFT. Neste caso comparar-se-á apenas o ganho de cada uma
das frequências de interesse. Foi preciso então calcular o erro entre os ganhos das frequências
do modelo e do sistema. Para isso foi criado um índice de erro entre os ganhos das frequências
do modelo e do sistema chamado ErroW . A equação (3.22) apresenta o cálculo deste índice.
3.4 Comentários Finais
47
O cálculo de ErroW trata-se de uma soma ponderada do erro de todas as frequências em
relação às frequências geradas pelo sistema elevadas ao quadrado. O cálculo realizado nessa
forma reforça a importância de representar prioritariamente as frequências de maior ganho e
secundariamente as frequências de menor ganho. Isso porque as frequências de maiores ganhos
tem maior participação na composição do sinal e consequentemente terão maior participação
na análise no domínio da frequência do sistema, enquanto às de menor ganho possuem menor
expressividade.
Logo, como a resposta em frequência diz respeito a dinâmica do sistema a ser modelado, isso
garante que o modelo obtido não seja apenas um ajuste aos dados de identif cação, mas, que
esteja incorporando a dinâmica presente no conjunto de dados.
3.4 Comentários Finais
Neste capítulo foi apresentada e discutida a metodologia desenvolvida e utilizada neste trabalho.
Foi discutido que o método apresentado pode tanto servir de método único de seleção de estrutura como servir de uma ferramenta de auxílio a métodos clássicos presentes na literatura
entre os quais pode-se citar o critério de ordenação ERR e o critério de informação de Akaike.
A grande vantagem de utilizar os métodos em conjunto seria diminuir o tempo de processamento e de convergência do algoritmo devido a grande redução do conjunto de agrupamentos
possíveis.
Nos dois próximos capítulos será apresentada a aplicação do método proposto em sistemas simulados e em sistemas reais. Serão avaliados o desempenho e a robustez do método submetido
à diferentes aplicações.
48
3 Seleção de estrutura através da análise de resposta em frequência
Capítulo 4
Estudo de casos simulados
4.1 Introdução
Neste capítulo a metodologia proposta neste trabalho é aplicada em dois sistemas simulados.
O primeiro é um sistema visto em (HABER; UNBEHAUEN, 1990). Este sistema foi escolhido principalmente por ser frequentemente utilizado na literatura para realização de identif cação. Isso permite que sejam comparados desempenhos de modelos obtidos anteriormente
com representações apontadas pelo método. Uma outra vantagem é poder obter sistemas com
não linearidades pré-estabelecidas variando um parâmetro do modelo.
O segundo sistema trata-se da simulação de uma planta de neutralização de pH. Seu uso neste
trabalho deve-se principalmente à sua característica não linear severa. Isso garante que o método
utilizado possa ser submetidos à casos extremos onde a identif cação nem sempre funciona utilizando ferramentas existentes na literatura. Casos recentes na literatura como em (CORRÊA,
2001; CAMPOS, 2007) apresentam estruturas de grande complexidade para uma representação satisfatória desse sistema. A metodologia aplicada fará uma busca por uma estrutura mais
compacta e que mesmo assim represente bem a dinâmica da planta de pH simulada.
4.2 Sistema simulado 1 - CTV
O sistema simulado 1 é um sistema de primeira ordem apresentado em (HABER; UNBEHAUEN, 1990):
50
4 Estudo de casos simulados
[1 + αu(t)]
dy
+ y(t) − u(t) = 0
dt
(4.1)
ou seja, trata-se de um sistema cujo ganho é constante e a constante de tempo varia com o ponto
de operação na seguinte forma τ (u) = 1 + αu(t). A equação (4.1) foi simulada com α = 5,
gerando conjuntos de dados distintos para identif cação e para validação assim como feito em
(CORRÊA; AGUIRRE, 2004). Não foi adicionado ruído para este caso.
Para a identif cação realizada no procedimento, foi usado um sinal com a entrada de amplitude
aleatória, excursionando o sistema em toda a sua faixa de operação com duração de cada patamar também aleatório. O tempo de amostragem usado foi de Ts = 0,1 s e a massa de dados
total possui 15000 observações.
Esse sistema foi simulado utilizando o Simulinkr em formato de diagrama de blocos utilizando
a equação (4.1). A simulação é baseada no algoritmo de Dormand-Prince chamado no programa
de ode45.m. A Figura 4.1 apresenta o diagrama de blocos da representação.
Figura 4.1: Sistema simulado 1 em representação por blocos.
Os dados gerados pela simulação podem ser visualizados nas Figuras 4.2 e 4.3, respectivamente.
Procurou-se excitar o sistema em toda a faixa de frequência onde o sistema tem potência espectral com o objetivo de excitar toda sua dinâmica.
4.2 Sistema simulado 1 - CTV
51
Dados de identificação: entrada vs. tempo
p.u.
4
3
2
1
0
0
100
0
100
200
300
400
500
600
tempo(s)
Dados de identificação: saída vs. tempo
200
300
700
3
p.u.
2.5
2
1.5
1
0.5
400
tempo(s)
500
600
700
Figura 4.2: Dados gerados pela simulação do sistema 1 para identif cação. Em (a), os dados de entrada
são apresentados e em (b) os dados de saída.
Dados de validação: entrada vs. tempo
p.u.
4
3
2
1
800
900
800
900
1000
1100
1200
1300
tempo(s)
Dados de identificação: saída vs. tempo
1400
3
p.u.
2.5
2
1.5
1
0.5
1000
1100
1200
tempo(s)
1300
1400
Figura 4.3: Dados gerados pela simulação do sistema 1 para validação. Em (a), os dados de entrada são
apresentados e em (b) os dados de saída.
É interessante notar que os vetores de dados de identif cação e de validação possuem tamanhos
iguais. Neste primeiro caso o sistema foi simulado, portanto, poderiam ser gerados quantos
52
4 Estudo de casos simulados
dados fossem necessários, não havendo a obrigatoriedade de critérios de divisão entre dados
identif cação e de validação. No entanto todos os critérios apreciados na seção 2.2.1.1 foram
observados. A frequência de amostragem utilizada na coleta foi de 62,83 rad/s.
Na Figura 4.4 a resposta em frequência do sistema, obtido de forma experimental, sendo o
módulo da razão entre a DTFT dos sinais de saída e entrada é apresentado.
Espectro de frequência do sistema simulado 1: CTV
10
0
Ganho (dB)
−10
−20
−30
−40
−50
−60
−2
10
−1
10
0
10
1
10
Frequencia (rad/s)
Figura 4.4: Ganho das frequências do sistema simulado 1 em dB. A linha reforça o comportamento do
ganho no domínio da frequências dos dados.
Pela Figura 4.4 pode-se observar que o comportamento do sistema simulado 1 é semelhante a
de um f ltro passa baixas e que, a frequência de corte desse f ltro considerando um ganho de
corte de -3dB, localiza-se entre 0,1 e 0,2 rad/s. Os pontos após esta faixa de frequência que
possuem ganhos que diferem do comportamento do sistema são causados pelo janelamento dos
dados. Foi utilizado a janela de Hammimg para diminuição desse efeito. A linha escura enfatiza
a relação entre os ganhos de entrada e saída do sistema simulado 1.
Diante das informações apresentadas nesta seção, pode-se partir para a obtenção da estrutura
seguindo a metodologia proposta.
4.2 Sistema simulado 1 - CTV
53
4.2.1 Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 1
utilizando o método clássico e a sugestão de grau de não linearidade proposta
Como já comentado, a frequência de corte do sistema simulado 1 encontra-se entre 0,1 e 0,2
rad/s. Logo foi escolhida a frequência 0,02 rad/s para excitar o sistema e verif car quantas e
quais são as frequências ressonantes.
A Figura 4.5 (a) apresenta parte do sinal gerado com frequência única de 0,02 rad/s com média
2,5 p.u. e amplitude máxima de 4,8 aplicado ao sistema simulado 1. A Figura 4.5 (b) apresenta
parte da resposta do sistema ao sinal aplicado.
(a) Parte dos dados de entrada para o sistema 1: CTV vs. tempo na frequência de 0,02 rad/s
p.u.
4
3
2
1
0
100
200
300
400
500
600
tempo(s)
(b) Parte dos dados do sistema 1: CTV saída vs. tempo submetido a uma frequência de 0,02 rad/s
0
100
200
p.u.
4
3
2
1
300
400
tempo(s)
500
600
Figura 4.5: Parte dos dados para o sistema 1 excitado a uma frequência de 0,02 rad/s. Em (a), os dados
de entrada são apresentados e em (b) os dados de saída.
A detecção visual de não linearidades a partir da análise dos sinais de entrada e saída é difícil
devido à semelhança deste conjunto de dados verif cado nas Figuras 4.5 (a) e (b) para entrada
e saída, respectivamente. Porém, o sinal parece ter sofrido uma pequena alteração na saída.
Como pode ser observado, houve uma mudança nos tempos de pico superior do sinal próximo
a 100 e a 400 segundos da saída em relação à entrada. Para uma análise mais criteriosa são
realizadas as DTFT’s tanto para o sinal de entrada quanto para o sinal de saída.
As Figuras 4.6 (a) e (b) apresentam a análise no domínio da frequência do sistema simulado 1
submetido a uma entrada de 0,02 rad/s. A análise no domínio da frequência foi obtida através
54
4 Estudo de casos simulados
da aplicação da DTFT tanto nos dados de entrada quanto nos dados de saída.
(a)
(b)
Espectro de frequência da entrada do sistema 1 submetido a uma frequência de 0,02 rad/s
Espectro de frequência da saída do sistema 1 submetido a uma frequência de 0,02 rad/s
2
2
1.5
|Y(jω)|
|U(jω)|
1.5
1
1
0.5
0.5
0
0
0.05
0.1
0.15
Frequencia (rad/s)
0.2
0.25
0
0
0.05
0.1
0.15
Frequencia (rad/s)
0.2
0.25
Figura 4.6: DTFT dos dados do sistema simulado 1 excitado a uma frequência de 0,02 rad/s. A DTFT do
sinal de entrada está representados em (a). A DTFT dos dados de saída estão representados
em (b). O eixo horizontal foi cortado na frequência 0,25 rad/s para uma melhor visualização.
Nas Figuras 4.6 pode-se verif car a existência de apenas uma frequência no sinal de entrada
e o surgimento de uma frequência ressonante no sinal de saída. A Tabela 4.1 apresenta as
frequências surgidas na saída e os seus ganhos em módulo.
Tabela 4.1: Frequências e respectivos ganhos surgidos na saída do sistema simulado 1.
Frequência (rad/s)
0,02
0,04
Ganho (módulo)
2,3911
0,1845
Analisando a Figura 4.6 (b) e confrontando com os dados da Tabela 4.1 observa-se o surgimento
de uma outra frequência aparente além da frequência de excitação. A primeira é a frequência
fundamental, resposta direta à entrada e principal componente da saída. Sugere-se então que
esta é a frequência que deve ser melhor representada pelo modelo. A única frequência ressonante do modelo é a de valor 0,04 rad/s. Esta frequência tem papel fundamental na dinâmica do
sistema 1 porque surge dentro de sua faixa de operação e dentro do limite de 30% da frequência
de amostragem.
Diante disso e da discussão feita na seção 3.3.1 sugere-se que o uso de um modelo polinomial de
grau de não-linearidade dois seria suf ciente para representação do sistema simulado 1. Graus de
4.2 Sistema simulado 1 - CTV
55
não linearidade maiores que este poderiam inserir frequências espúrias do sistema no modelo.
Os agrupamentos de termos que podem ser sugeridos para a estrutura do modelo, analisando
pela resposta em frequência, são: os termos lineares, capazes de reproduzir a frequência fundamental dos dados da saída e em conjunto representar o comportamento de f ltro passa-baixas do
sistema; o agrupamento de termos cruzados de grau de não linearidade par. Estes últimos agrupamentos colaboram para a representação da frequência ressonante múltipla par da frequência
fundamental; e um agrupamento de termo puro de entrada de grau de não linearidade 2. Este
colabora para a reprodução da segunda frequência ressonante exibida pelo sistema.
Os máximos atrasos na entrada e saída foram determinados utilizando o FCC.
O modelo obtido pelo método clássico da literatura é apresentado na equação (4.2).
y(k) = 1,6680y(k − 1) − 0,6795y(k − 2) + 0,0119u(k − 1) − 0,0008u(k − 1)2
+0,3294u(k − 1)y(k − 2) − 0,3287u(k − 1)y(k − 1) − 0,0005u(k − 2)
+6,42 × 10−5 u(k − 2)2 − 2,08 × 10−5 u(k − 2)u(k − 1) − 1,3343y(k − 1)2
−1,3154y(k − 2)2 + 2,6499y(k − 2)y(k − 1) − 0,0005
(4.2)
O modelo da equação (4.2) apresentou 13 termos em sua estrutura ordenados do mais relevante
para o menos relevante na explicação da dinâmica do sistema. O RMSE foi de 0,1290. A
predição livre e a análise no domínio da frequência do modelo são apresentadas nas Figuras 4.7
(a) e (b) respectivamente.
56
4 Estudo de casos simulados
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
3.5
2
1.8
3
1.6
2.5
1.4
1.2
p.u.
|Y(jω)|
2
1
1.5
0.8
0.6
1
0.4
0.5
0
700
Dados
Modelo (4.2)
800
900
0.2
1000
1100
1200
tempo(s)
1300
1400
1500
0
0
0.02
0.04
0.06
Frequencia (rad/s)
0.08
0.1
Figura 4.7: Predição livre e análise no domínio da frequência do modelo obtido pelo método clássico
para o sistema simulado 1. A predição livre está ilustrada em (a). A análise no domínio da
frequência em (b).
Na Figura 4.7 pode-se observar que a análise no domínio da frequência obtida pelo modelo
apresentam as frequências fundamental e outras duas ressonantes. Os valores das frequências
ressonantes para este caso são apresentados na Tabela 4.2.
Tabela 4.2: Frequências e respectivos ganhos surgidos na saída modelo obtido pelo método clássico para
o sistema simulado 1.
Frequência (rad/s)
0,02
0,04
0,06
Ganho (módulo)
2,076
0,2248
0,0045
Pela Tabela 4.2 é possível observar que as frequências apresentadas pela predição livre do sistema submetido ao sinal de frequência única de 0,02 rad/s são as mesmas apresentadas pelo modelo exceto a segunda frequência ressonante. A frequência 0,06 rad/s apresentada pelo modelo
não corresponde ao conjunto de frequências ressonantes apresentadas pelo sistema. Portanto, o
surgimento desta frequência pode representar uma dinâmica espúria afetando o desempenho do
modelo. O ErroW calculado foi de 0,1233.
Diante disso, procurou-se aplicar o método proposto para verif car a possibilidade de melhorar
o desempenho do modelo e ainda avaliar a necessidade de cada termo e agrupamento de termos
apresentado na estrutura fornecida pelo método clássico.
4.2 Sistema simulado 1 - CTV
57
4.2.2 Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 1
utilizando a resposta em frequência
Nesta seção será implementada a metodologia proposta no capítulo 3 com o objetivo de melhorar a representação apresentada pelo método clássico.
Analisando os termos apresentados pelo modelo da equação (4.2) em relação à resposta em
frequência, é possível detectar agrupamentos que contribuem de forma semelhante para a resposta em frequência do modelo. Um exemplo disso, seria a contribuição para o surgimento da
primeira ressonante apresentada pelo sistema. De acordo com a Tabela 3.1 os agrupamentos de
termos puros de entrada com grau de não linearidade par, os agrupamentos de termos puros de
saída de grau de não linearidade par e os agrupamentos de termos cruzados contribuem para o
surgimento desta segunda frequência ressonante no modelo.
Portanto, alguns desses termos poderiam ser descartados do modelo fornecido pelo método
clássico por representarem a mesma contribuição na resposta em frequência do sistema. Considerando a ordem de importância apresentada pelo critério ERR os termos e agrupamentos que
devem ser considerados são: y(k − 1), y(k − 2), u(k − 1) para representação da dinâmica
principal do modelo como um f ltro passa-baixas e representação da frequência fundamental.
E também, o termo u(k − 1)2 , responsável pela representação da frequência ressonante apresentada pelo sistema. Com esta estrutura, estima-se novamente os parâmetros e obtém-se o
seguinte modelo:
y(k) = 1,8075y(k − 1) − 0,8099y(k − 2) + 0,0027u(k − 1) − 0,0003u(k − 1)2 (4.3)
A Figura 4.8 apresenta em (a) a predição livre do modelo obtido e em (b) a análise no domínio
da frequência. O RMSE foi de 0,9904.
58
4 Estudo de casos simulados
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
3.5
1
3
2.5
0.8
p.u.
|Y(jω)|
2
0.6
1.5
0.4
1
0.2
0.5
0
700
Dados
Modelo (4.3)
800
900
1000
1100
1200
tempo(s)
1300
1400
1500
0
0
0.02
0.04
0.06
Frequencia (rad/s)
0.08
0.1
Figura 4.8: Predição livre e análise no domínio da frequência do primeiro modelo obtido para o sistema
simulado 1. A predição livre está representados em (a). A análise no domínio da frequência
em (b).
A Tabela 4.3 apresenta os valores dos ganhos das frequências fornecidas pelo modelo apresentado pelo método no passo 1.
Tabela 4.3: Frequências e respectivos ganhos surgidos na saída modelo obtido pelo primeiro passo do
método proposto para o sistema simulado 1.
Frequência (rad/s)
0,02
0,04
Ganho (módulo)
1,1552
0,3322
Pode-se observar através da Figura 4.8 e pelo índice RMSE que o desempenho do modelo foi
inferior ao do modelo obtido pelo método clássico. Na Tabela 4.3 nota-se que os ganhos de cada
frequência fornecida pelo modelo têm valores bem diferentes daqueles obtidos pelo sistema
apresentados na Tabela 4.1. O ErroW calculado foi de 0,4833 bem maior que o ErroW
apresentado pelo modelo apresentado na equação (4.2).
Diante do pior desempenho do modelo da equação (4.3) comparado ao modelo obtido pelo
método clássico, sugere-se acrescentar mais um termo que possa colaborar nos ganhos das
frequências fundamental e ressonante exibidas na saída do modelo. Seguindo a ordenação realizada pelo critério ERR, adiciona-se o termo u(k − 1)y(k − 2) porque este termo também
possui contribuição para a frequência ressonante apresentada pelo sistema. Depois estima-se
novamente os parâmetros e realiza-se uma nova predição livre. O modelo obtido é apresentado
4.2 Sistema simulado 1 - CTV
59
na equação (4.4).
y(k) = 1,6923y(k − 1) − 0,6995y(k − 2) + 0,0064u(k − 1) − 0,0012u(k − 1)2
(4.4)
+0,0014u(k − 1)y(k − 2)
A Figura 4.9 apresenta em (a) a predição livre do modelo obtido e em (b) a análise no domínio
da frequência . O RMSE foi de 0,3690.
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
3.5
2
3
2.5
1.5
p.u.
|Y(jω)|
2
1.5
1
1
0.5
0.5
0
700
Dados
Modelo (4.4)
800
900
1000
1100
1200
tempo(s)
1300
1400
1500
0
0
0.02
0.04
0.06
Frequencia (rad/s)
0.08
0.1
Figura 4.9: Predição livre e e análise no domínio da frequência do segundo modelo obtido para o sistema
simulado 1. A predição livre está representados em (a). A análise no domínio da frequência
em (b).
A Tabela 4.4 apresenta os valores dos ganhos das frequências fornecidas pelo modelo apresentado pelo método no passo 2.
Tabela 4.4: Frequências e respectivos ganhos surgidos na saída modelo obtido pelo segundo passo do
método proposto para o sistema simulado 1.
Frequência (rad/s)
0,02
0,04
0,06
0,08
Ganho (módulo)
2,205
0,2856
0,119
0,04821
60
4 Estudo de casos simulados
Pode-se notar pelas informações da Figura 4.1 (b) e da Tabela 4.4 que pelo efeito da adição
de um termo cruzado ao modelo surgem outras frequências na saída fornecida pelo modelo.
Embora isso tenha ocorrido, o modelo teve seu desempenho signif cativamente melhorado se
comparado ao modelo (4.2). Isso pode ter sido causado pela aproximação dos ganhos das
frequências de interesse (0,02 e 0,04 rad/s) do modelo aos ganhos das frequências reais do
sitema. O ErroW obtido foi de 0,4833.
No próximo passo, utilizando o critério de ordenação ERR, o agrupamento que será acrescentado é u(k − 1)y(k − 1). Com este agrupamento o modelo obtido e estimado é apresentado na
equação (4.5).
y(k) = 1,7851y(k − 1) − 0,7951y(k − 2) + 0,0099u(k − 1) − 0,0004u(k − 1)2
(4.5)
+0,3795u(k − 1)y(k − 2) − 0,3791u(k − 1)y(k − 1)
A Figura 4.10 apresenta em (a) a pedição livre do modelo obtido e em (b) a análise no domínio
da frequência . O RMSE foi de 0,1293.
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
3.5
2
1.8
3
1.6
2.5
1.4
1.2
p.u.
|Y(jω)|
2
1.5
1
0.8
0.6
1
0.4
0.5
0
700
Dados
Modelo (4.5)
800
900
1000
0.2
1100
1200
tempo(s)
1300
1400
1500
0
0
0.02
0.04
0.06
Frequencia (rad/s)
0.08
0.1
Figura 4.10: Predição livre e análise no domínio da frequência do terceiro modelo obtido para o sistema
simulado 1. A predição livre está representados em (a). A análise no domínio da frequência
em (b).
A Tabela 4.5 apresenta os valores dos ganhos das frequências fornecidas pelo modelo apresentado pelo método no passo 3.
4.2 Sistema simulado 1 - CTV
61
Tabela 4.5: Frequências e respectivos ganhos surgidos na saída modelo obtido pelo terceiro passo do
método proposto para o sistema simulado 1.
Frequência (rad/s)
0,02
0,04
0,06
0,08
Ganho (módulo)
2,019
0,2123
0,038
0,009
Comparando os valores dos ganhos das frequências 0,02 e 0,04 rad/s das Tabelas 4.1 e 4.5
é possível notar que estes continuam próximos. Isso signif ca que as frequências que devem
ser representadas pelo modelo tem o seu ganho próximos aos ganhos reais do sistema. No
entanto, em relação às frequências espúrias ao sistema 1, 0,06 e 0,08 rad/s, a adição do termo
u(k − 1)y(k − 1) causou uma forte atenuação em seus ganhos fazendo com que a resposta
em frequência do modelo fosse melhor aproximada da resposta em frequência do sistema. O
ErroW calculado foi de 0,1449 ainda maior que o modelo (4.5).
No próximo passo, utilizando o critério de ordenação ERR, o termo que será acrescentado é
u(k − 2). Este é um termo linear e contribui apenas para o ganho da frequência fundamental. Tanto pode ser descartado como pode ser utilizado para um novo cálculo de ErroW .
Isso porque a sua contribuição pode ser irrelevante já que o ganho da frequência fundamental
encontra-se muito próximo ao ganho da frequência do sistema. Também pode ter efeito danoso
na representação da dinâmica principal do sistema. Com o propósito de análise, será apresentada a predição livre também para o modelo com este termo. Com este termo o modelo obtido
e estimado é apresentado na equação (4.6).
y(k) = 1,7846y(k − 1) − 0,7946y(k − 2) + 0,0101u(k − 1) − 0,0004u(k − 1)2
+0,3793u(k − 1)y(k − 2) − 0,3789u(k − 1)y(k − 1) − 0,0003u(k − 2) (4.6)
A Figura 4.11 apresenta em (a) a predição livre do modelo obtido e em (b) a análise no domínio
da frequência. O RMSE foi de 0,1301.
62
4 Estudo de casos simulados
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
3.5
2
1.8
3
1.6
2.5
1.4
1.2
p.u.
|Y(jω)|
2
1.5
1
0.8
0.6
1
0.4
0.5
0
700
Dados
Modelo (4.6)
800
900
1000
0.2
1100
1200
tempo(s)
1300
1400
1500
0
0
0.02
0.04
0.06
Frequencia (rad/s)
0.08
0.1
Figura 4.11: Predição livre e análise no domínio da frequência da saída do quarto modelo obtido para o
sistema simulado 1. A predição livre está representados em (a). A análise no domínio da
frequência em (b).
A Tabela 4.6 apresenta os valores dos ganhos das frequências fornecidas pelo modelo apresentado pelo método no passo 4.
Tabela 4.6: Frequências e respectivos ganhos surgidos na saída modelo obtido pelo quarto passo do
método proposto para o sistema simulado 1.
Frequência (rad/s)
0,02
0,04
0,06
0,08
Ganho (módulo)
2,02
0,2125
0,038
0,009
Confrontando os dados das Tabelas 4.5 e 4.6 pode-se notar que a adição do termo u(k − 2)
não teve inf uência signif cativa na resposta em frequência do modelo. Por isso, poderia ser
dispensado da análise no método proposto. o ErroW obtido para este caso foi de 0,1445.
Descartando o termo u(k − 2) e adicionando o termo u(k − 2)2 é realizado o próximo passo. O
modelo obtido agora é apresentado na equação (4.7).
4.2 Sistema simulado 1 - CTV
63
y(k) = 1,7849y(k − 1) − 0,7949y(k − 2) + 0,0098u(k − 1) − 0,0004u(k − 1)2
+0,3794u(k − 1)y(k − 2) − 0,3790u(k − 1)y(k − 1)
(4.7)
−4,24 × 10−5 u(k − 2)2
A Figura 4.12 apresenta em (a) a predição livre do modelo obtido e em (b) a sua análise no
domínio da frequência. O RMSE foi de 0,1301.
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
3.5
2
1.8
3
1.6
2.5
1.4
1.2
p.u.
|Y(jω)|
2
1.5
1
0.8
0.6
1
0.4
0.5
0
700
Dados
Modelo (4.7)
800
900
1000
0.2
1100
1200
tempo(s)
1300
1400
1500
0
0
0.02
0.04
0.06
Frequencia (rad/s)
0.08
0.1
Figura 4.12: Predição livre e análise no domínio da frequência do quinto modelo obtido para o sistema
simulado 1. A predição livre está representados em (a). A DTFT em (b).
A Tabela 4.7 apresenta os valores dos ganhos das frequências fornecidas pelo modelo apresentado pelo método no passo 5.
Tabela 4.7: Frequências e respectivos ganhos surgidos na saída modelo obtido pelo quinto passo do
método proposto para o sistema simulado 1.
Frequência (rad/s)
0,02
0,04
0,06
0,08
Ganho (módulo)
2,019
0,2124
0,038
0,009
Novamente para este quinto modelo a alteração da resposta em frequência ocasionada pela
adição do termo u(k − 2)2 foi pequena. O ErroW calculado para a DTFT deste quinto modelo
64
4 Estudo de casos simulados
foi de 0,1448. Pode-se notar que a adição dos dois últimos termos não alterou signif cativamente
a resposta em frequência da saída do modelo. O cálculo de ErroW praticamente não se alterou.
A Figura 4.13 apresenta a evolução de ErroW ao longo dos passos.
Evolução de ErroW
0.5
0.45
0.4
ErroW
0.35
0.3
0.25
0.2
0.15
0.1
1
1.5
2
2.5
3
3.5
4
4.5
Passo (Número de agrupamentos de termos adicionados)
5
Figura 4.13: Evolução do cálculo de ErroW para cada passo para o sistema simulado 1.
Nota-se que nos passos quatro e cinco não há diminuição signif cativa do ErroW calculado.
Isso sugere que tanto o termo u(k − 2) e u(k − 2)2 são irrelevantes na explicação da dinâmica
do sistema no domínio da frequência. Portanto pode-se haver o corte da estrutura no passo três.
O modelo sugerido pelo método será o modelo apresentado na equação (4.5).
É importante salientar que essa avaliação não é feita para os primeiros termos lineares. Isso
porque esses termos não colaboram sozinhos na reprodução da frequência ressonante da saída
e não alcançam desempenho satisfatório na representação de sistemas que apresentam não linearidades.
Agrupamento de termos que geram as mesmas frequências podem ser sugeridos ao modelo. A
contribuição de cada agrupamento para cada frequência dependerá também do coef ciente que
acompanha o agrupamento. Para o termo u(k − 1)2 é fácil predizer o efeito do termo sobre o
modelo. Porém para os termos cruzados, u(k − 1)y(k − 2) e u(k − 1)y(k − 1) esta conclusão
não é tão trivial, pois, estes termos podem fornecer inf nitas frequências ressonantes quando
adicionados ao modelo. No entanto, diferentemente de u(k − 1)2 , os agrupamentos cruzados contribuem de forma signif cativa para a frequência fundamental do sistema representado.
4.3 Sistema simulado 2 - A planta de neutralização de pH simulada
65
Porém, o grande problema em adicionar estes agrupamentos, é o fornecimento de frequências
espúrias à resposta em frequência do sistema. Quando o termo u(k − 1)y(k − 2) é adicionado
ao modelo há um melhora razoável em seu desempenho, porém pela análise da resposta em frequência, surgem frequências não pertencentes à dinâmica do sistema. A adição do agrupamento
u(k − 1)y(k − 1) parece ter efeito contrário nas frequências ressonantes funcionando como um
f ltro que diminui os seus ganhos.
Uma outra análise que pode ser feita é quanto à importância do agrupamento u(k − 1)2 sugerida
pelo ERR. Entre os termos não lineares, o ERR o classif ca como o mais signif cativo na representação da dinâmica do sistema. Isso corrobora a sugestão desse agrupamento para compor o
modelo pela análise da resposta em frequência da saída do sistema, pois, é o termo que contribuirá diretamente para o ganho da primeira frequência ressonante do sistema em discussão.
Comparando o modelo obtido através do método clássico e o modelo obtido pelo método proposto neste trabalho pode-se notar que o RMSE calculado é bastante semelhante. Para o modelo
clássico o erro RMSE foi de 0,1290 e para o modelo proposto foi de 0,1293. Como esses erros
são praticamente iguais pode-se af rmar que o desempenho dos modelos são equivalentes. Só
que há uma vantagem do modelo obtido pelo método proposto em relação ao modelo obtido
pelo método clássico. A relevante redução da complexidade do modelo. Vários agrupamentos de termos foram descartados diminuindo também a possibilidade de o modelo representar
dinâmicas espúrias ao sistema.
4.3 Sistema simulado 2 - A planta de neutralização de pH
simulada
O segundo sistema simulado trata-se de uma planta de neutralização de pH implementado por
(CAMPOS, 2007). A simulação foi feita com base em um modelo fenomenológico baseado na
aproximação físico-química de invariantes de reação. A escolha desse sistema se dá pelo fato de
apresentar uma complexa e difícil não linearidade que pode variar sensivelmente sob pequenas
mudanças nas condições do processo.
Simplif cadamente a planta de neutralização de pH simulada têm três entradas e uma saída:
• Fluxo de entrada de um ácido forte;
• Fluxo de entrada de uma base forte;
• Fluxo de entrada de uma solução tampão;
66
4 Estudo de casos simulados
• Medição de pH.
Neste trabalho a entrada única trata-se da vazão de reagente básico e a saída, também única,
trata-se da medição de pH. As demais entradas são mantidas constantes. Melhores informações
sobre a planta de neutralização de pH simulada podem ser obtidas em (CAMPOS, 2007).
Para a identif cação realizada no procedimento, foi usado um sinal com a entrada de amplitude
aleatória, excursionando o sistema em toda a sua faixa de operação com duração de cada patamar também aleatório. O tempo de amostragem usado foi de Ts = 15 s e a massa de dados total
possui 48000 observações.
Parte dos dados de identif cação e validação podem ser visualizados nas Figuras 4.14 e 4.15
respectivamente.
Dados de identificação: entrada vs. tempo
4
p.u.
3
2
1
0
1
2
3
4
5
6
tempo(s)
Dados de identificação: saída vs. tempo
7
4
x 10
10
p.u.
8
6
4
0
1
2
3
4
tempo(s)
5
6
7
4
x 10
Figura 4.14: Parte dos dados gerados pela simulação do sistema 2 para identif cação. Em (a), os dados
de entrada são apresentados e em (b) os dados de saída.
4.3 Sistema simulado 2 - A planta de neutralização de pH simulada
67
Dados de validação: entrada vs. tempo
4
p.u.
3
2
1
3.6
3.7
3.8
3.9
4
4.1
4.2
tempo(s)
Dados de identificação: saída vs. tempo
4.3
5
x 10
10
p.u.
8
6
4
3.6
3.7
3.8
3.9
4
tempo(s)
4.1
4.2
4.3
5
x 10
Figura 4.15: Parte dos dados gerados pela simulação do sistema 2 para validação. Em (a), os dados de
entrada são apresentados e em (b) os dados de saída.
Espectro de frequência do sistema simulado 2: pH
30
20
Ganho (dB)
10
0
−10
−20
−30
−40
−4
10
−3
10
Frequencia (rad/s)
−2
10
Figura 4.16: Ganho das frequências do sistema simulado 2 em dB. A linha escura reforça o comportamento do ganho no domínio da frequência dos dados.
Para a simulação com uma frequência def nida foram aplicados sinais senoidais pré-def nidos
com amplitude f xa também excursionando toda a faixa de interesse do sistema com o objetivo
68
4 Estudo de casos simulados
de varrer toda sua dinâmica e consequentemente excitar todas as suas não linearidades. Na
Figura 4.16 a resposta em frequência obtida com base nos dados de identif cação são apresentados.
Assim como no sistema simulado 1, o comportamento do sistema simulado 2 aproxima-se ao
comportamento de um f ltro passa-baixas. O corte de frequências, considerando um ganho de
-3dB localiza-se entre 4,0 × 10−3 e 6,0 × 10−3 rad/s. A linha escura reforça a indicação do
comportamento no domínio da frequência do sistema simulado 2.
Obtidas as informações descritas até aqui, parte-se para a seleção da estrutura do modelo para
o sistema simulado 2.
4.3.1 Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 2
utilizando o método clássico e a sugestão de grau de não linearidade proposta
Uma vez conhecida a região de frequência de corte do sistema simulado 2 pode-se estabelecer
a frequência que será usada para excitá-lo. A frequência escolhida foi de 2,0 × 10−4 rad/s.
Tal frequência encontra-se dentro da região de potência espectral do sistema 2 e é cerca de dez
vezes inferior a mais alta frequência da região de potência espectral do sistema.
A Figura 4.17 (a) apresenta o sinal de entrada utilizado na excitação do sistema 1. Este sinal
possui apenas uma componente de frequência. A média estabelecida é de 2,5 e amplitude
máxima de 5. Já a Figura 4.17 (b) apresenta um sinal de saída diferente da senoidal aplicada.
Isso conf rma a presença de frequências diferentes no sinal de saída em relação ao sinal de
entrada.
4.3 Sistema simulado 2 - A planta de neutralização de pH simulada
69
(a)
Parte dos dados de entrada para o sistema 2 vs. tempo na frequência de 2,0× 10−3 rad/s
5
p.u.
4
3
2
1
0
2
3
4
5
6
7
4
tempo(s)
x 10
(b)
Parte dos dados do sistema 2 saída vs. tempo submetido a uma frequência de 2,0× 10−3 rad/s
2
3
4
10
p.u.
8
6
4
5
6
7
tempo(s)
4
x 10
Figura 4.17: Parte dos dados para o sistema 2 excitado a uma frequência de 2,0 × 10−4 rad/s. Em (a), os
dados de entrada são apresentados e em (b) os dados de saída.
A Figura 4.18 (a) e (b) apresenta a análise no domínio da frequência do sistema simulado 2
submetido a uma entrada de 2,0 × 10−4 rad/s. A a análise no domínio da frequência foi obtida
através da aplicação da DTFT tanto nos dados de entrada quanto nos dados de saída.
(a)
(b)
Espectro de frequência da saída (y) para o sistema 2: pH
−3
submetido a uma frequência de 2,0× 10 rad/s
Espectro de frequência da entrada (u) para o sistema 2: pH
−3
submetido a uma frequência de 2,0× 10 rad/s
1.2
2
1
0.8
|Y(jω)|
|U(jω)|
1.5
0.6
1
0.4
0.5
0.2
0
0
0
0.5
1
1.5
2
2.5
3
Frequencia (rad/s)
3.5
4
4.5
5
−3
x 10
0
0.5
1
1.5
2
2.5
3
Frequencia (rad/s)
3.5
4
4.5
5
−3
x 10
Figura 4.18: Análise no domínio da frequência dos dados do sistema simulado 2 excitado a uma frequência 2,0 × 10−4 rad/s. A análise no domínio da frequência do sinal de entrada está
representados em (a). Para os dados de saída, em (b).
70
4 Estudo de casos simulados
Nas Figuras 4.18 pode-se comprovar a existência de apenas uma frequência da entrada e o
surgimento de outras frequências na saída. Diferentemente do sistema simulado 1, o espectro
de frequência do sinal de saída do sistema simulado 2 apresenta diversas frequências na saída.
Isso ocorre, provavelmente, em decorrência da severa não-linearidade da planta de neutralização
de pH simulada.
O desaf o neste ponto é determinar quais frequências são relevantes na formação do sinal de
saída quando o sistema é submetido a um sinal de frequência única igual a 2,0 × 10−4 rad/s e
quais frequências que o modelo será capaz de reproduzir em sua saída. Sugere-se pelo desenvolvimento realizado que um grau de não-linearidade extremamente elevado possa representar
bem o sistema no domínio da frequência. Porém, isto é pouco interessante na prática pois torna
o modelo muito complexo fazendo com que sua simulação e implementação f que trabalhosa.
Mesmo que não se tenha um grau de não linearidade próximo do número de frequências ressonantes do sistema, pela Tabela 3.1 pode-se observar que agrupamentos de termos puros de
saída e agrupamentos de termos cruzados fornecem ao sistema inf nitas frequências ressonantes.
Logo, existe a possibilidade de se adotar graus de não linearidade menores, diminuindo a complexidade do modelo e, ao mesmo tempo, possibilitando o surgimento de inf nitas frequências
ressonantes no modelo acrescentando os agrupamentos mencionados. Percebe-se também que
modelos com apenas agrupamentos puros de entrada não serão capazes de reproduzir bem a
dinâmica do sistema.
Através dessa discussão e da discussão realizada no capítulo 3, foi escolhido um grau de não
linearidade três para obter o modelo pelo método clássico.1 Sabe-se que o grau de não linearidade não atende à reprodução de todas as frequências na saída, portanto, é certo que serão
necessários agrupamento de termos cruzados e talvez agrupamentos de termos puros de saída.
Foram determinados máximos atrasos de entrada e saída um e três respectivamente através do
FCC.
Uma outra questão é em relação a quantidade de frequências que deverão ser consideradas no
cálculo do índice ErroW . Para tanto, obteve-se o modelo pelo método clássico e observou-se a
sua capacidade de representação da resposta em frequência. Através dessa análise determinouse quais frequências seriam utilizadas para determinar o cálculo de ErroW . Contudo, para este
caso, o método procurará uma estrutura que apresente um melhor desempenho sendo menos
complexa e partirá do modelo obtido pelo método clássico.
O modelo obtido pelo método clássico é apresentado na equação (4.8).
1
Testes suprimidos neste trabalho mostraram que o aumento do grau de não linearidade para quatro ou cinco
não melhoram satisfatoriamente o desempenho do modelo para os dados do sistema simulado 2.
4.3 Sistema simulado 2 - A planta de neutralização de pH simulada
71
y(k) = 1,0424y(k − 1) − 0,6490y(k − 2) + 0,0772u(k − 1) + 0,4129y(k − 3)
−0,0497u(k − 1)2 y(k − 1) + 0,3492u(k − 1)y(k − 1)
−0,0279u(k − 1)y(k − 3)2 + 0,2965 + 0,0342u(k − 1)3
+0,0312u(k − 1)2 y(k − 3) − 0,2261u(k − 1)y(k − 3) − 0,1404u(k − 1)2
+0,0530u(k − 1)y(k − 3)y(k − 1) − 0,0271u(k − 1)y(k − 1)2
(4.8)
+0,003y(k − 3)2
A ordem dos termos apresentada na equação (4.8) é a ordem hierárquica fornecida pelo ERR.
O critério de informação de Akaike forneceu um valor de 32 termos para o modelo, porém o
mesmo não convergiu. O ponto de corte foi def nido em uma busca exaustiva retirando-se os
termos até encontrar uma estrutura para que o modelo convergisse. Chegou-se ao número de
termos igual a 15.
O RMSE do modelo apresentado na equação (4.8) foi de 0,2099. A predição livre e a resposta
em frequência do modelo obtido pelo ERR para o sistema simulado 2 são apresentadas nas
Figuras 4.19 (a) e (b), respectivamente.
(a)
(b)
Parte dos dados de validação: saída
Espectro de frequência da saída (y)
4.5
10
4
9
3.5
8
3
7
pH
|Y(jω)|
6
5
4
2.5
2
1.5
3
1
2
Dados
Modelo (4.8)
1
0
3.64
3.66
3.68
3.7
3.72 3.74
tempo(s)
0.5
3.76
3.78
3.8
3.82
5
x 10
0
0
0.5
1
1.5
Frequencia (rad/s)
2
2.5
−3
x 10
Figura 4.19: Predição livre e análise no domínio da frequência do modelo obtido pelo método clássico
para o sistema simulado 2. A predição livre está representada em (a). A análise no domínio
da frequência em (b).
A Figura 4.19 (b) apresenta todas as frequências que o modelo foi capaz de criar na saída
submetido ao mesmo sinal de entrada de frequência única aplicado ao sistema. A Tabela 4.8
72
4 Estudo de casos simulados
apresenta as frequências surgidas na saída fornecida pelo modelo apresentado na equação (4.8)
submetido ao sinal de entrada de frequência única igual a 2,0 × 10−4 rad/s.
Tabela 4.8: Frequências e respectivos ganhos surgidos na saída do modelo obtido pelo método clássico
para o sistema simulado 2.
Frequência (rad/s)
2,0 × 10−4
4,0 × 10−4
6,0 × 10−4
8,0 × 10−4
1,0 × 10−3
1,2 × 10−3
1,4 × 10−3
Ganho (módulo)
4,511
0,5131
0,4681
0,2601
0,2045
0,1276
0,0796
A Tabela 4.8 mostra que o modelo obtido pelo ERR foi capaz de representar a frequência
fundamental e mais seis frequências ressonantes com ganhos consideráveis. Portanto, para
o cálculo de ErroW , serão consideradas as primeiras sete frequências surgidas na saída do
sistema para realização do corte de termos. A Tabela 4.9 apresenta as sete primeiras frequências
e seus respectivos ganhos para o sistema simulado 2 observadas na Figura 4.18 (b).
Tabela 4.9: Frequências e respectivos ganhos surgidos na saída do sistema simulado 2.
Frequência (rad/s)
2,0 × 10−4
4,0 × 10−4
6,0 × 10−4
8,0 × 10−4
1,0 × 10−3
1,2 × 10−3
1,4 × 10−3
Ganho (módulo)
2,475
0,206
0,497
0,1145
0,1469
0,0423
0,0319
Com base nas Tabelas 4.8 e 4.9 pode-se calcular o ErroW para a o modelo obtido pelo ERR.
Para este, o valor calculado de ErroW é 0,5884.
4.3 Sistema simulado 2 - A planta de neutralização de pH simulada
73
4.3.2 Obtenção de estrutura de modelo polinomial NARX para o sistema simulado 2
utilizando a resposta em frequência
Diferentemente da forma apresentada para o sistema 1, a obtenção do modelo pelo método
proposto neste capítulo será simplif cada para o sistema simulado 2. Os modelos não serão
simulados a cada retirada de termos ou de agrupamentos. Apenas será apresentada na Tabela
4.10 os valores de ErroW para cada um dos passos apresentados. Os valores de ErroW para
os termos lineares localizados nas primeiras posições classif cados pelo ERR são considerados
como número de agrupamentos igual a zero na Tabela.
O objetivo do método neste segundo caso será o de melhorar a estrutura do modelo obtido pelo
ERR com 15 termos de forma que alcance um melhor desempenho ou diminua sua complexidade mantendo o desempenho. Com isso, o método proposto funcionará como um critério
de informação, ou de corte, retirando agrupamentos e termos que não explicam a dinâmica do
sistema, ou ainda, incorporam dinâmicas espúrias ao modelo.
Tabela 4.10: Valores de ErroW para cada adição de termo no modelo para o sistema simulado 2.
N.o Ω adicionado ErroW
0
1,0225
1
0,8943
2
0,8126
3
0,7852
4
0,7912
5
0,6655
6
0,5044
7
0,4972
8
0,4955
0,6050
9
Na Figura 4.20 é apresentado o ErroW a cada agrupamento acrescentado ao modelo. Nota-se
que a partir do passo 7 no qual é acrescentado o termo u(k − 1)y(k − 3) não há grande variação
do índice ErroW . Portanto, pode-se dizer que o acréscimo de termos a partir de então não
tem efeito na resposta em frequência do sistema. Portanto, limita-se a estrutura do modelo no
sétimo agrupamento após os termos lineares.
74
4 Estudo de casos simulados
Evolução de ErroW
1.1
1
ErroW
0.9
0.8
0.7
0.6
0.5
0.4
1
2
3
4
5
6
7
8
9
Passo (Número de agrupamentos de termos adicionados)
10
Figura 4.20: Evolução do cálculo de ErroW para cada passo para o sistema simulado 2.
O modelo sugerido pelo método é apresentado na equação (4.9).
y(k) = 1,1824y(k − 1) − 0,6404y(k − 2) − 0,2164u(k − 1) + 0,3180y(k − 3)
−0,0630u(k − 1)2 y(k − 1) + 0,2989u(k − 1)y(k − 1)
−0,0014u(k − 1)y(k − 3)2 + 0,3703 + 0,0180u(k − 1)3
+0,0449u(k − 1)2 y(k − 3) − 0,1921u(k − 1)y(k − 3)
(4.9)
O modelo 4.9 apresenta 11 termos enquanto o modelo obtido pela busca exaustiva foi estimado
com 15 termos. É importante salientar que o número de termos fornecido pelo critério de corte
de Akaike forneceu um número de termos de um modelo que não convergiu (31 termos). A
predição livre do modelo 4.9 é apresentada na Figura 4.21 (a) e em (b) a sua análise no domínio
da frequência submetido ao sinal de frequência única de 2,0 × 10−4 rad/s.
4.3 Sistema simulado 2 - A planta de neutralização de pH simulada
(a)
75
(b)
Parte dos dados de validação: saída
Espectro de frequência da saída (y)
4
10
3.5
9
8
3
7
2.5
|Y(jω)|
pH
6
5
4
2
1.5
3
1
2
Dados
Modelo (4.9)
1
0
3.64
3.66
3.68
3.7
3.72 3.74
tempo(s)
0.5
3.76
3.78
3.8
3.82
5
x 10
0
0
0.5
1
1.5
Frequencia (rad/s)
2
2.5
−3
x 10
Figura 4.21: Predição livre e análise no domínio da frequência do modelo obtido pelo método de análise
da resposta em frequência para o sistema simulado 2. A predição livre está representados
em (a). A análise no domínio da frequência em (b).
O método proposto apresentou um modelo que obteve um desempenho satisfatório. O índice
RMSE calculado foi de 0,2099, o mesmo índice calculado para o modelo com 15 termos. Isso
sugere que os quatro termos a mais presentes no modelo apresentado na equação (4.8) não
interferiam na representação da dinâmica dos dados apresentados.
O cálculo de ErroW não é necessário para os termos lineares hierarquizados pelo ERR. Isso
porque esses termos só colaboram para a reprodução da frequência fundamental não registrando
interferência nas dinâmicas não lineares do sistema. E neste trabalho, procurou-se trabalhar com
os agrupamentos de termos que explicam as dinâmicas não lineares dos sistemas representados.
Embora, os índices de ErroW para modelos lineares, provavelmente, teriam grandes módulos
por apresentarem ganho zero nas frequências ressonantes.
4.3.3 Validação estática dos modelos
Como discutido na seção 2.2.1.7 uma das formas de validar modelos estimados é através da
curva estática. Muito embora a curva estática não explicite a dinâmica do sistema, seu comportamento deve aproximar bem da curva estática do sistema modelado para que todos os seus
pontos f xos possam ser aproximados pelo modelo.
Na Figura 4.22 são apresentados os comportamentos estáticos exibidos pelo sistema e pelos
modelos obtidos neste trabalho para a planta simulada de neutralização de pH.
76
4 Estudo de casos simulados
Curva Estática − planta (pH)
11
10
9
pH
8
7
6
Dados
5
Modelo (4.9) − Mét. proposto
Modelo (4.8) − Mét. clássico
4
3
2
0.5
1
1.5
2
2.5
Vazão de base (ml/s)
3
3.5
4
Figura 4.22: Curva estática do sistema simulado 2 e dos modelos obtidos pelo critério de ERR e pelo
método proposto. A linha contínua apresenta a curva estática do sistema simulado 2. A linha
descontínua mais clara apresenta a característica estática do modelo obtido pelo ERR com
15 termos. A linha descontínua mais escura apresenta a característica estática do modelo
obtido pelo método proposto baseado na resposta em frequência.
Na Figura 4.22 pode-se observar que a curva estática do sistema simulado 2 é melhor aproximada quando se segue a sugestão do método proposto. Isso pode explicar que os termos que
foram retirados do modelo não têm importância signif cativa na representação do sistema. Além
disso, sugere-se que o modelo com os termos a mais, embora tenha tido desempenho semelhante quando é comparado o índice RMSE, fornece um comportamento estático bem diferente
do comportamento estático real do sistema. Enquanto o modelo apresentado pelo método proposto é capaz de apresentar nesta mesma região um comportamento muito mais próximo do
sistema. Isso pode ser observado nas regiões de vazão de base entre 0 e 0,5 ml/s.
No entanto, como pode ser visto na f gura 4.14 poucos dados de identif cação estão presentes
nas regiões de melhora de desempenho do comportamento estático do modelo (4.9). Portanto,
não se pode af rmar que a melhora nesse desempenho é referente à aplicação do método.
4.4 Comentários Finais
77
4.4 Comentários Finais
Neste capítulo foram expostos e discutidos os resultados da aplicação do método em dois casos simulados. O primeiro caso, trata-se de um sistema com uma não linearidade mais suave
enquanto no segundo caso a não linearidade é mais severa. Embora o sistema simulado 2 ter
caráter não-linear mais severo, observou-se que o método obteve um considerável desempenho
na seleção de estrutura do modelo NARX.
A utilização do método proposto como critério de corte sugeriu modelos menos complexos para
os dois casos e ainda assim, manteve-se um desempenho satisfatório.
No próximo capítulo o método proposto será usado na seleção de estrutura de modelos para
sistemas reais.
78
4 Estudo de casos simulados
Capítulo 5
Estudo de casos reais
5.1 Introdução
No Capítulo 4 foi apresentada a aplicação do método proposto no Capítulo 3 em dois sistemas
simulados. Neste capítulo a metodologia é aplicada em sistemas reais e por isso apresentam
algumas peculiaridades. Entre essas características peculiares dos sistemas reais, pode-se citar
o aparecimento de ruído nos sinais coletados e a dif culdade em obter dados em grandes quantidades, seja pelo tempo requerido para simulação ou pelas condições não ideais de operação,
que podem não permitir a obtenção de dados que varram toda a faixa de operação do sistema,
prejudicando assim a identif cação.
Neste capítulo, o método proposto no Capítulo 3 foi implementado em dois sistemas reais. O
primeiro sistema trata-se de um sistema térmico de bancada em que sua fonte de calor é uma
lâmpada incandescente. O segundo sistema trata-se de um modelo de um helicóptero de dois
graus de liberdade.
A validação dos modelos foi feita pelo método tradicional utilizando dados separados para a
validação diferente do conjunto de dados de identif cação e foi escolhido o sistema térmico para
se realizar a validação do modelo obtido pelo comportamento estático.
5.2 Sistema Térmico
A principal característica de um sistema térmico é a conservação de calor em seu interior. Isso
faz com que o aumento de temperatura medida no interior do sistema térmico seja mais rápido
80
5 Estudo de casos reais
enquanto a diminuição de temperatura, ao cessar a fonte de calor, seja mais lenta. Devido à esta
variação de constante de tempo, o sistema térmico apresenta um caráter não linear.
Resumidamente, o sistema térmico usado neste trabalho trata-se de uma montagem de bancada
composta por um tubo de PVC de 200mm de diâmetro e 300mm de altura apoiado em uma superfície de madeira e uma lâmpada incandescente em seu interior. O forno é fracamente isolado
termicamente de forma que variações da temperatura ambiente afetam o seu comportamento
dinâmico. Os dados consistem em uma entrada de controle variando de 5 a 0 Volts que são
eletronicamente convertidos em tensão aplicada à uma lâmpada incandescente. Esta lâmpada é
de 60 W de potência e é alimentada por uma tensão variável de 0 a 127 Volts. Quanto maior a
tensão aplicada à lâmpada, maior é a energia térmica fornecida ao sistema. A saída do sistema
é exibida em uma escala de 0 a 10 Volts e é fornecida por um circuito eletrônico que contém
como sensor um transistor LM35. Cada variação unitária de tensão de saída representa a variação de cerca de 10o C na temperatura do sistema. Os dados são trocados entre o sistema e o
computador através de uma placa USB-6008 da National Instrumentsr
Nas seções subsequentes serão apresentadas as formas de aquisição de dados e a aplicação do
método no auxílio à obtenção de uma estrutura adequada de modelos polinomiais.
5.2.1 Aquisição de dados
Os dados de identif cação devem ser capazes de ref etir os diversos regimes dinâmicos apresentados pelo sistema. Para tanto, o sinal de excitação da planta deve apresentar espectro suf cientemente amplo em frequência e amplitude de tal forma que excursione o sistema pelos regimes
dinâmicos de interesse. Infelizmente, isso nem sempre é possível em sistemas reais, devido às
restrições operacionais.
Os dados foram coletados a uma taxa de amostragem de 62,83 rad/s. A f m de evitar mal
condicionamento numérico na estimação de parâmetros devido a medições repetitivas, os dados
foram dizimados para uma taxa de amostragem de 6,2832 rad/s para realização do procedimento
de identif cação.
Embora para a identif cação, dados repetitivos possam implicar em mal condicionamento numérico,
na etapa de estimação de parâmetros, para a obtenção da análise no domínio da frequência,
optou-se por manter a taxa de amostragem inicial dos dados a f m de melhorar a visualização
das frequências de interesse.
As Figuras 5.1 e 5.2 apresentam os dados coletados e separados em dados para a identif cação
e dados para validação.
5.2 Sistema Térmico
81
Tensão (Volts)
(a) Dados de identificação: entrada vs. tempo
4
3
2
1
1
1.5
2
2.5
3
3.5
tempo(s)
(b) Dados de identificação: saída vs. tempo
4
4
x 10
Temperatura (Cº)
70
60
50
40
1
1.5
2
2.5
tempo(s)
3
3.5
4
4
x 10
Figura 5.1: Dados gerados pelo sistema térmico para identif cação. Em (a), os dados de entrada são
apresentados e em (b) os dados de saída.
Tensão (Volts)
(a) Dados de validação: entrada vs. tempo
4
3
2
Temperatura (Cº)
1
0
1000
0
1000
2000
3000
4000
5000
6000
tempo(s)
(b) Dados de identificação: saída vs. tempo
60
55
50
45
40
2000
3000
4000
tempo(s)
5000
6000
Figura 5.2: Dados gerados pelo sistema térmico para validação. Em (a), os dados de entrada são apresentados e em (b) os dados de saída.
Na Figura 5.3 a resposta em frequência do sistema é apresentado.
82
5 Estudo de casos reais
Espectro de frequência do sistema real 1: Forno
40
30
20
Ganho (dB)
10
0
−10
−20
−30
−40
−3
10
−2
−1
10
10
0
10
Frequencia (rad/s)
Figura 5.3: Ganho das frequências do sistema real térmico em dB. A linha escura reforça o comportamento do ganho no domínio da frequências dos dados.
A resposta em frequência apresentada pela Figura 5.3 demonstra que o sistema térmico comportase como um f ltro passa-baixas atenuando fortemente frequências acima da faixa aproximada
de 2,0 × 10−2 e 6,0 × 10−2 rad/s. Fisicamente, um sistema térmico excitado a um sinal de alta
frequência não é capaz de responder ao sinal de entrada. Isso porque o processo tem característica dinâmica lenta. Portanto, sinais de alta frequência inseridos na entrada são fortemente
atenuados na saída.
A linha escura reforça a visualização do resposta em frequência do sistema térmico mais aproximada.
5.2.2 Obtenção de estrutura de modelo polinomial NARX para o sistema térmico utilizando o método clássico e a sugestão de grau de não linearidade proposta
A frequência de corte de operação do sistema térmico encontra-se próxima a faixa de 2,0×10−2
e 6,0 × 10−2 rad/s como pode ser visto na Figura 5.3. Conhecida a região de espectro de
frequência onde o sistema térmico funciona, determinou-se a frequência de excitação do sistema
para verif cação do surgimento de frequências ressonantes em 5,3 × 10−3 rad/s.
A Figura 5.4 (a) apresenta o sinal de entrada aplicado ao sistema térmico e a resposta do sistema
ao sinal aplicado na Figura 5.4 (b) à frequência determinada.
5.2 Sistema Térmico
83
Temperatura (Cº)
(a)
Parte dos dados de entrada para o forno vs. tempo na frequência de 0,0053 rad/s
4
3
2
1
1
1.05
1.1
1.15
1.2
tempo(s)
1.25
1.2
tempo(s)
1.25
1.3
1.35
4
x 10
(b)
Parte dos dados do forno saída vs. tempo submetido a uma frequência de 0,0053 rad/s
Temperatura (Cº)
65
60
55
50
45
1
1.05
1.1
1.15
1.3
1.35
4
x 10
Figura 5.4: Parte dos dados para o sistema térmico excitado a uma frequência de 5,3 × 10−3 rad/s. Em
(a), os dados de entrada são apresentados e em (b) os dados de saída.
(a)
(b)
Espectro de frequência da entrada (u) para o forno submetido a uma frequência de 0,0053 rad/s
Espectro de frequência da saída (y) para o forno submetido a uma frequência de 0,0053 rad/s
9
1.6
8
1.4
7
1.2
6
|Y(jω)|
|U(jω)|
1
0.8
0.6
5
4
3
0.4
2
0.2
1
0
0
0.05
0.1
0.15
Frequencia (rad/s)
0.2
0.25
0
0
0.05
0.1
0.15
Frequencia (rad/s)
0.2
0.25
Figura 5.5: Análise no domínio da frequência do sistema térmico excitado a uma frequência de 5,3×10−3
rad/s. A análise no domínio da frequência do sinal de entrada está representados em (a) e dos
dados de saída estão representados em (b).
As Figuras 5.5 (a) e (b) apresentam a análise no domínio da frequência do sistema térmico
submetido a uma entrada de 5,3 × 10−3 rad/s da entrada e da saída respectivamente. A análise
84
5 Estudo de casos reais
no domínio da frequência foi obtida através da aplicação da DTFT tanto nos dados de entrada
quanto nos dados de saída.
Nas Figuras 5.5 (a) e (b) pode-se comprovar a existência de apenas uma frequência da entrada
e o surgimento de duas frequências ressonantes existentes na saída, muito embora haja o aparecimento de frequências devido à amostragem e ao janelamento do sinal. A Tabela 5.1 apresenta
as frequências ressonantes da frequência fundamental surgidas na saída e os seus ganhos em
módulo.
Tabela 5.1: Frequências e respectivos ganhos surgidos na saída do sistema térmico.
Frequência (rad/s)
5,3 × 10−3
1,06 × 10−2
1,59 × 10−2
Ganho
8,577
0,4084
0,5522
Os dados exibidos na Tabela 5.1 indicam uma não linearidade moderada nos dados do sistema
térmico devido ao ganho das ressonâncias serem baixos em relação a frequência aplicada na
entrada. No entanto, para representar de forma adequada o sistema como sugere o método,
pode-se utilizar grau de não-linearidade dois ou três. Optou-se pelo menor grau de não linearidade a f m de obter um modelo menos complexo. Logo, agrupamentos de termos cruzados
ou puros de entrada deverão estar presentes para representar a terceira frequência ressonante já
que agrupamentos puros de entrada de grau de não linearidade par não as reproduzem. Quanto
aos atrasos máximos de entrada e saída foram determinados, pelo FCC, dois para os dois casos.
Por ser um sistema real e estar sujeito a ruídos, obteve-se modelos NARMAX com 20 termos
de ruído lineares.
O modelo obtido pelo método clássico é apresentado na equação (5.1)
y(k) = 3,5972y(k − 1) − 2,5138y(k − 2) − 0,2710u(k − 1)y(k − 1) − 0,3879y(k − 1)2
−0,3541y(k − 2)2 + 0,0048u(k − 2)y(k − 2) + 0,2561u(k − 1)y(k − 2)
+0,7412y(k − 2)y(k − 1) − 0,0290u(k − 2)y(k − 2) + 0,4167u(k − 1)
20
X
−1,5933 +
θ̂i ξ(k − 1) + ξ(k)
(5.1)
i=1
O modelo da equação (5.1) apresentou 11 termos em sua estrutura ordenados do mais relevante
5.2 Sistema Térmico
85
para o menos relevante na explicação da dinâmica do sistema e os termos de ruído. O rmse foi
de 0,3401. A predição livre e a resposta em frequência do modelo são apresentadas nas Figuras
5.6 (a) e (b), respectivamente.
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
14
62
60
12
58
10
54
|Y(jω)|
Temperatura (Cº)
56
52
50
8
6
48
46
4
Dados
Modelo (5.1)
44
2
42
1
1.1
1.2
1.3
1.4
tempo(s)
1.5
0
1.6
0
0.01
4
x 10
0.02
0.03
Frequencia (rad/s)
0.04
0.05
Figura 5.6: Predição livre e análise no domínio da frequência do modelo obtido pelo método clássico
para o sistema térmico. A predição livre está representados em (a). A análise no domínio da
frequência em (b).
Pela Figura 5.6 pode-se observar o surgimento de frequências devido ao janelamento do sinal e
da amostragem. No caso de sistemas reais esse efeito é intensif cado pois os conjuntos de dados
têm tamanhos limitados. Um outro problema que surge é a presença de ruído nos dados. Esses
ruídos aparecem nas respostas em frequência do sistema e assim como o ajanelamento do sinal,
podem camuf ar as frequências ressonantes que surgem na saída com menores ganhos causando
erros no cálculo do índice ErroW para implementação do método.
A Tabela 5.2 apresenta as frequências e seus respectivos ganhos para o sistema térmico observados na Figura 5.6 (b).
Tabela 5.2: Frequências e respectivos ganhos surgidos na saída do modelo obtido pelo método clássico
para o sistema térmico.
Frequência (rad/s)
5,3 × 10−3
1,06 × 10−2
1,59 × 10−2
Neste caso, o valor de ErroW foi 0,6123.
Ganho (módulo)
14,39
0,5117
0
86
5 Estudo de casos reais
5.2.3 Obtenção de estrutura de modelo NARMAX polinomial para o sistema térmico
utilizando a resposta em frequência
Tendo ordenado os termos do modelo hierarquicamente utilizando o ERR, aplica-se o método
proposto para a realização do corte da estrutura. Os valores de ErroW são apresentados na
Tabela 5.6.
Tabela 5.3: Valores de ErroW para cada adição e termo no modelo para o sistema térmico.
N.o Ω adicionado ErroW
1
0,8856
2
0,6165
3
0,6325
0,5305
4
5
0,5634
6
0,6269
7
0,7276
Na Figura 5.7 é apresentado o ErroW a cada agrupamento acrescentado ao modelo. Note que
a partir do passo 4 no qual é acrescentado o termo u(k − 2)y(k − 2) o índice ErroW volta a
crescer. Portanto, pode-se dizer que o acréscimo de termos a partir de então tem efeito negativo
na resposta em frequência do sistema, ou seja, a inclusão de termos está incluindo dinâmicas
espúrias ao modelo. Portanto, limita-se a estrutura do modelo no quarto agrupamento após os
termos lineares.
Evolução de ErroW
0.9
0.85
0.8
ErroW
0.75
0.7
0.65
0.6
0.55
0.5
1
2
3
4
5
6
Passo (Número de agrupamentos de termos adicionados)
7
Figura 5.7: Evolução do cálculo de ErroW para cada passo para o sistema térmico.
5.2 Sistema Térmico
87
O modelo sugerido pelo método é apresentado na equação (5.2).
y(k) = 2,0014y(k − 1) − 0,9726y(k − 2) − 0,0088u(k − 1)y(k − 1) − 0,0060y(k − 1)2
20
X
2
+0,0056y(k − 2) + 0,0056u(k − 2)y(k − 2) +
θ̂i ξ(k − 1) + ξ(k)
(5.2)
i=1
Observe que o modelo (5.2) tem apenas seis termos, cinco a menos que o modelo (5.1), apresentado pelo método clássico. A predição livre do modelo obtido é apresentada na Figura
5.8.
Dados de validação: saída
62
60
58
Temperatura (Cº)
56
54
52
50
48
46
Dados
Modelo (5.2)
44
42
1
1.1
1.2
1.3
1.4
tempo(s)
1.5
1.6
4
x 10
Figura 5.8: Predição livre do modelo obtido pelo método proposto para o sistema térmico.
O método proposto apresentou o modelo 5.2 que obteve um desempenho satisfatório. O índice
RMSE calculado foi de 0,3560 próximo ao RMSE obtido pelo modelo (5.1).
5.2.4 Validação estática dos modelos
O sistema térmico foi escolhido para a realização da validação estática.
Na Figura 5.9 são apresentados os comportamentos estáticos exibidos pelo sistema e pelos
modelos obtidos neste trabalho para o sistema térmico.
88
5 Estudo de casos reais
Curva Estática − forno
Modelo ERR+AKAIKE
Modelo Resp. em frequência
80
Temperatura (Cº)
70
60
50
40
30
0.5
1
1.5
2
Tensão (Volts)
2.5
3
Figura 5.9: Curva estática do sistema térmico e dos modelos obtidos pelo critério de ERR+AKAIKE
e pelo método proposto. A linha contínua apresenta a curva estática do sistema térmico.
A linha descontínua mais clara apresenta a característica estática do modelo obtido pelo
ERR+AKAIKE. A linha mais escura apresenta a característica estática do modelo obtido
pelo método proposto baseado na resposta em frequência.
Pela Figura pode-se observar que as curvas estáticas dos dois modelos aproximaram de forma
insuf ciente à curva estática do sistema térmico. No entanto, quando é comparado o modelo
5.2 pelo método proposto e o modelo (5.1), obtido pelo método clássico, observa-se que os
modelos têm desempenho estático próximos. Embora não tenha tido um desempenho superior
signif cativo, o método proposto sugere um modelo mais compacto o que o torna vantajoso
perante o outro modelo.
5.3 O helicóptero de 2 graus de liberdade
O sistema real 2 usado neste trabalho trata-se de um protótipo simplif cado de um helicóptero
fabricado pela empresa Quanserr . Este sistema é considerado uma aproximação de um helicóptero real porque possui apenas 2 GDL (Graus de Liberdade): a arfagem e a guinada. A
Figura 5.10 apresenta os eixos de rotação do helicóptero para o movimento de arfagem e de
guinada. A arfagem constitui o movimento de rotação em torno do eixo horizontal e a guinada
é o movimento de rotação em torno do eixo vertical.
5.3 O helicóptero de 2 graus de liberdade
89
Figura 5.10: Os dois eixos de rotação do helicóptero.
As duas entradas deste sistema são independentes e constituem-se em tensão de controle de 0
a 24 Volts convertidas em tensão de alimentação do motor de cada uma das hélices. Tanto a
localizada acima do helicóptero quanto àquela localizada na cauda conforme a foto do protótipo
apresentada na Figura 5.11.
Figura 5.11: Foto do protótipo de helicóptero com dois GDL.
As duas saídas do protótipo constituem-se nos ângulos de cada um dos movimentos de rotação.
A primeira, o ângulo de rotação do movimento de arfagem e a segunda o ângulo de rotação para
o movimento de guinada.
Este trabalho limitou-se ao uso de sistemas SISO. Portanto a entrada do motor da cauda foi
mantida constante e a leitura do ângulo de rotação do movimento de guinada foi desprezada.
Mais detalhes sobre a coleta de dados serão apresentados na seção 5.3.1.
A ação da gravidade contribui para uma condição especial ao movimento do sistema. Ao mesmo
tempo que faz oposição ao movimento de subida, colabora para o movimento de descida inf uindo na velocidade de subida e descida do movimento de arfagem. Isso faz com que o
sistema tenha ao menos duas constantes de tempo. A partir disso já pode-se considerar que o
90
5 Estudo de casos reais
sistema apresenta caráter não-linear.
5.3.1 Aquisição de dados
Como escrito anteriormente, a leitura do ângulo de rotação do movimento de guinada foi desprezada mantendo-se apenas a leitura do ângulo de rotação do movimento de arfagem. A entrada
desse sistema trata-se da alimentação do motor da hélice superior.
Os dados foram coletados a uma taxa de amostragem de 62,83 rad/s. A f m de evitar mal
condicionamento numérico na estimação de parâmetros devido a medições repetitivas, os dados
foram decimados para uma taxa de amostragem de 6,2832 rad/s para realização do procedimento de identif cação. As Figuras 5.12 e 5.13 apresentam os dados usados na identif cação e
os separados para a validação, respectivamente.
Muito embora os dados de entrada do sistema real 2 permitissem uma aplicação de tensão na
faixa de 0 a 24 Volts, optou-se por usar a faixa de 13 a 15 Volts na entrada do sistema para
que o helicóptero não atingisse os limites de movimento superior ou inferior do movimento de
rotação de guinada.
5.3.2 Obtenção de estrutura de modelo NARMAX polinomial para o helicóptero utilizando o método clássico e a sugestão de grau de não linearidade proposta
A frequência de 0,0628 rad/s foi escolhida para compor o sinal de frequência única aplicado ao
sistema e aos modelos na implementação do método.
A Figura 5.14 (a) apresenta o sinal de entrada aplicado ao protótipo e a resposta do sistema ao
sinal aplicado na Figura 5.14 (b) à frequência determinada.
5.3 O helicóptero de 2 graus de liberdade
91
Tensão (Volts)
(a) Dados de identificação: entrada vs. tempo
14.5
14
13.5
Ângulo − Arfagem (º)
13
0
200
400
600
800
tempo(s)
(b) Dados de identificação: saída vs. tempo
0
200
400
1000
0
−10
−20
−30
−40
600
tempo(s)
800
1000
Figura 5.12: Dados gerados pelo helicóptero para identif cação. Em (a), os dados de entrada são apresentados e em (b) os dados de saída.
Tensão (Volts)
(a) Dados de validação: entrada vs. tempo
14.5
14
13.5
Ângulo − Arfagem (º)
1100
1150
1200
1250
1300
1350
tempo(s)
(b) Dados de identificação: saída vs. tempo
−10
−15
−20
−25
−30
−35
1100
1150
1200
1250
tempo(s)
1300
1350
Figura 5.13: Dados gerados pelo helicóptero para validação. Em (a), os dados de entrada são apresentados e em (b) os dados de saída.
92
5 Estudo de casos reais
Tensão (Volts)
(a) Dados de entrada para o helicóptero vs. tempo na frequência de 0.6283 rad/s
14.5
14
13.5
400
500
600
700
800
tempo(s)
(b) Dados do helicóptero saída vs. tempo submetido a uma frequência de 0.6283 rad/s
Ângulo − Guinada (º)
100
200
300
100
200
300
−10
−15
−20
−25
−30
400
500
tempo(s)
600
700
800
Figura 5.14: Parte dos dados para o helicóptero excitado a uma frequência de 0,0628 rad/s. Em (a), os
dados de entrada são apresentados e em (b) os dados de saída.
A Figura 5.5(a) e (b) apresenta a análise no domínio da frequência do sinal de arfagem do
helicóptero submetido a uma entrada de 0,0628 rad/s. A resposta em frequência foi obtida
através da aplicação da FFT tanto nos dados de entrada quanto nos dados de saída.
(a)
(b)
Espectro de frequência do sinal de arfagem do helicóptero
submetido a uma frequência de 0.6283 rad/s
Espectro de frequência do sinal de arfagem do helicóptero
submetido a uma frequência de 0.6283 rad/s
8
3
7
2.5
6
2
|Y(jω)|
|U(jω)|
5
1.5
4
3
1
2
0.5
1
0
0
0.1
0.2
0.3
0.4
Frequencia (rad/s)
0.5
0.6
0.7
0
0
0.1
0.2
0.3
0.4
Frequencia (rad/s)
0.5
0.6
0.7
Figura 5.15: Análise no domínio da frequência dos dados do sinal de arfagem do helicóptero excitado a
uma frequência de 0,0628 rad/s. A análise no domínio da frequência dos dados de entrada
está representados em (a) e dos dados de saída estão representados em (b).
5.3 O helicóptero de 2 graus de liberdade
93
Nas Figuras 5.15 (a) e (b) a visualização das frequências está mais complexa que nos sistemas
anteriores. Pode-se observar que a inf uência do janelamento, mesmo minimizado com os f ltros
de janela de Hamming, dos dados de entrada insere diversas frequências na DTFT da entrada.
Não obstante, essas frequências acabam por interferir na resposta em frequência do sistema
obtida pela DTFT do sinal de saída. Porém, é possível visualizar em (a) que a frequência
predominante é aquela aplicada na entrada de 0,0628 rad/s. Em (b), é possível ver a frequência
ressonante 0,1257 rad/s com ganho menor em relação à frequência fundamental. As demais
frequências surgidas podem ser em decorrência da discretização do sinal e do número limitado
de dados coletados. A Tabela 5.4 apresenta as frequências surgidas na saída e os seus ganhos
em módulo.
Tabela 5.4: Frequências e respectivos ganhos surgidos na saída do helicóptero.
Frequência (rad/s)
0,0628
0,1257
Ganho
8,545
0,835
Os dados exibidos na Tabela 5.4 indicam uma não linearidade moderada nos dados do helicóptero devido ao ganho das frequências ressonantes serem baixos em relação a frequência
aplicada na entrada. Para representar de forma adequada o sistema como sugere o método,
pode-se utilizar grau de não-linearidade dois pois há apenas mais uma frequência ressonante
visível. Quanto aos atrasos máximos de entrada e saída foram determinados em dois pelo FCC.
Por ser um sistema real e estar sujeito a ruídos, foram obtidos modelos NARMAX com 20
termos de ruído lineares.
O modelo obtido pelo método clássico é apresentado na equação (5.3).
y(k) = 1,9980y(k − 1) − 1,0004y(k − 2) − 0,0001y(k − 1)2 − 0,0020u(k − 2)y(k − 1)
−3,0961 + 0,0122u(k − 2)2 + 0,0008u(k − 1)2 + 0,0025y(k − 2)y(k − 1)
20
X
2
−0,0015y(k − 2) + 0,0021u(k − 2)y(k − 2) +
θ̂i ξ(k − 1) + ξ(k)
(5.3)
i=1
O modelo da equação (5.3) apresentou 10 termos em sua estrutura ordenados do mais relevante
para o menos relevante na explicação da dinâmica do sistema. O RMSE foi de 0,6843. A
predição livre e a análise no domínio da frequência do modelo são apresentadas nas Figuras
5.16 (a) e (b) respectivamente.
94
5 Estudo de casos reais
(a)
(b)
Dados de validação: saída
Espectro de frequência da saída (y)
8
Dados
Modelo (5.3)
−12
−14
7
−16
5
−20
|Y(jω)|
Ângulo − Arfagem (º)
6
−18
−22
4
−24
3
−26
−28
2
−30
1
−32
1145
1150
1155
1160
1165
tempo(s)
1170
1175
1180
0
0
0.02
0.04
0.06
Frequencia (Hz)
0.08
0.1
Figura 5.16: Predição livre e análise no domínio da frequência do modelo obtido pelo método clássico
para o helicóptero. A predição livre está representados em (a). A análise no domínio da
frequência em (b).
Pela Figura 5.16 pode-se observar o surgimento de frequências devido ao janelamento do sinal e
da amostragem. Esse efeito é intensif cado ainda mais pois os conjuntos de dados têm tamanhos
limitados. Os ruídos presentes no conjunto de dados aparecem nas respostas em frequência do
sistema e podem camuf ar as frequências ressonantes que surgem na saída com menores ganhos
dif cultando o cálculo de ErroW para implementação do método.
A Tabela 5.5 apresenta as frequências e seus respectivos ganhos para o helicóptero observadas
na Figura 5.6 (b). Neste caso, o valor de ErroW foi 0,0635.
Tabela 5.5: Frequências e respectivos ganhos surgidos na saída do modelo obtido pelo método clássico
para o helicóptero.
Frequência (rad/s)
0,0628
0,1257
Ganho
9,019
1,196
5.3.3 Obtenção de estrutura de modelo polinomial NARX para o helicóptero utilizando
a resposta em frequência
Na seção anterior, a análise no domínio da frequência do sinal de arfagem do protótipo do
helicóptero e um modelo pelo método clássico com 10 termos em sua estrutura foram obtidos.
Nesta seção, o método proposto será aplicado para efetuar o corte na estrutura do modelo. A
5.3 O helicóptero de 2 graus de liberdade
95
Tabela 5.6 mostra os valores de ErroW calculados a cada passo realizado pelo método.
Tabela 5.6: Valores de ErroW para cada adição e termo no modelo para o sistema térmico.
N.o Ω adicionado ErroW
1
0,9153
2
0,0900
3
0,0690
4
0,0816
5
0,0772
6
0,0750
7
0,0635
Na Figura 5.17 é apresentado o ErroW a cada agrupamento acrescentado ao modelo. Nota-se
que a partir do passo 3 no qual é acrescentado o termo u(k − 2)2 , o índice ErroW volta a
crescer. Portanto, pode-se dizer que o acréscimo de termos a partir de então tem efeito negativo
na resposta em frequência do sistema, ou seja, a inclusão de termos está incluindo dinâmicas
espúrias ao modelo. Portanto limita-se a estrutura do modelo no terceiro agrupamento após os
termos lineares.
Evolução de ErroW
1
0.9
0.8
0.7
ErroW
0.6
0.5
0.4
0.3
0.2
0.1
0
1
2
3
4
5
6
Passo (Número de agrupamentos de termos adicionados)
7
Figura 5.17: Evolução do cálculo de ErroW para cada passo para o helicóptero.
O modelo sugerido pelo método proposto é apresentado na equação (5.4).
96
5 Estudo de casos reais
y(k) = 2,1862y(k − 1) − 0,9475y(k − 2) + 0,0018y(k − 1)2 − 0,0145u(k − 2)y(k − 1)
20
X
−0,0836 +
θ̂i ξ(k − 1) + ξ(k)
(5.4)
i=1
Observe que o modelo (5.4) tem apenas cinco termos, cinco a menos que o modelo apresentado
pelo método clássico no modelo (5.3). A predição livre do modelo obtido é apresentada na
Figura 5.18.
Dados de validação: saída
Dados
Modelo (5.4)
−12
−14
Ângulo − Arfagem (º)
−16
−18
−20
−22
−24
−26
−28
−30
−32
1145
1150
1155
1160
1165
tempo(s)
1170
1175
1180
Figura 5.18: Predição livre do modelo obtido pelo método proposto para o helicóptero.
O método proposto apresentou o modelo (5.4) que obteve um desempenho fraco porém parecido
com o desempenho do modelo (5.3), obtido pelo método clássico. O índice RMSE calculado
foi de 0,7125 enquanto para o método clássico foi de 0,6843.
Pode-se observar que o índice de erro de simulação foi maior para o modelo obtido pelo
método mostrando um desempenho fraco em relação ao modelo obtido pelo método clássico.
A dif culdade em encontrar um modelo com uma estrutura que apresentasse um melhor desempenho pode ter sido causada pela interferência de frequências espúrias ao sistema surgidas
pelo ajanelamento do sinal no cálculo da FFT.
5.4 Comentários Finais
97
5.4 Comentários Finais
Dois sistemas reais foram utilizados neste capítulo para validar a ef ciência do método proposto.
Um sistema térmico foi escolhido pela sua facilidade de implementação e observação e pela
característica não linear apresentada assim, como o helicóptero, aquisição recente do grupo de
pesquisa MOCP.
A principal dif culdade apresentada na implementação da metodologia proposta foi a limitação
na quantidade de dados e a presença de ruídos. Essas ocorrências causaram interferências consideráveis na obtenção das DTFT’s de cada um dos sistemas, prejudicando assim, a análise
correta das respostas em frequência.
Resultados foram apresentados e comentados. Observou-se, contudo, que assim como em sistemas simulados, o método teve razoável desempenho fornecendo estruturas mais compactas e
ainda fornecendo indícios de termos de importância considerável na representação de cada um
dos sistemas.
98
5 Estudo de casos reais
Capítulo 6
Considerações Finais
A etapa de seleção de estrutura em identif cação de sistemas ainda não foi completamente solucionada. Diversos trabalhos têm sido realizados no intuito de facilitar o processo de busca de
uma estrutura mais próxima do ideal de um modelo sem a necessidade da experiência do modelador e se possível realizá-lo em posse de apenas o conjunto de dados em alguns casos. É
importante ressaltar que cada um tem a sua limitação.
Um critério amplamente discutido na literatura é o ERR. Apesar de sua ampla aplicabilidade,
tal critério pode escolher termos incorretos ou redundantes em condições não ideais de identif cação. Essas condições não ideais de identif cação referem-se a dados superamostrados ou
ruidosos ou quando o sinal de entrada é relativamente lento.
A f m de investigar o problema de seleção de estrutura de modelos NARX polinomiais, foi
proposto neste trabalho um estudo de comparação entre a análise da resposta em frequência de
um sistema e cada um dos agrupamentos de termos pertencentes a um modelo que o represente.
Uma análise matemática foi feita para cada tipo de agrupamento e algumas simulações foram
realizadas comprovando a aplicabilidade da metodologia.
No Capítulo 3 foi apresentada, analisada e discutida as relações de cada um dos termos possíveis em um determinado modelo. Foi visto que apenas para agrupamentos de termos puros
de entrada pode-se determinar quais as frequências geradas ao se aplicar um sinal de entrada
de frequência única no modelo. Também foi visto que agrupamento de termos cruzados ou
agrupamentos puros de saída fornecem inf nitas frequências ao modelo quando submetido ao
mesmo sinal. Pode-se perceber que este comportamento permitia a diminuição dos graus de
não linearidade de modelos pois para a representação de frequências ressonantes bastava adicionar um agrupamento de termos puro de saídas ou cruzados para gerar inf nitas ressonantes.
100
6 Considerações Finais
Observou-se que isso poderia justif car ser desnecessário graus de não linearidades muito altos
em alguns casos.
A partir dessas análises foi possível sugerir a inclusão ou exclusão de alguns termos de um modelo mesmo já hierarquizado pelo ERR, pois tornou-se possível avaliar o efeito do agrupamento
na resposta em frequência do modelo.
Os testes realizados para os sistemas simulados mostraram que a metodologia proposta funcionou como critério de corte de estruturas com agrupamentos hierarquizados pelo ERR. Observou-se que, em comparação ao critério de informação de Akaike, o desempenho dos modelos
sugeridos pelo método proposto foram semelhantes, porém, com a vantagem de fornecer estruturas mais compactas. Isso demonstrou que em alguns casos o critério de Akaike pode fornecer
termos que não explicam a dinâmica dos dados, algumas vezes fornecendo estruturas extremamente complexas. O método proposto foi capaz de reconhecer estruturas bem compactas que
tiveram desempenhos próximos a um modelo com uma estrutura muito mais complexa. A aplicação do método na planta simulada de pH mostrou que, embora no erro de simulação, não
houvesse melhora no desempenho dinâmico do modelo, ainda assim com uma estrutura bem
mais compacta, houve uma melhora signif cativa na aproximação do comportamento estático.
Para os casos reais, o método proposto não apresentou robustez à presença de ruído e ao ajanelamento do sinal. Isso porque o cálculo da FFT realizado para uma quantidade limitada de dados
causa o aparecimento de frequências não pertencentes ao sistema que camuf am as frequências
verdadeiras. A presença de ruídos nos dados também faz surgir frequências que dependendo
de sua amplitude esconde frequências que compõem o espectro real do sistema. Porém, modelos obtidos pelo método tiveram desempenhos semelhantes à modelos obtidos pelo método
clássico. Isso comprova que o método clássico também está sujeito à falhas na obtenção de
modelos devido à condições não ideais de identif cação.
O método aqui implementado serve como um critério de informação, porém, não é fundamentalmente estatístico como o critério de Akaike, mas apresenta relação com a representação da
dinâmica do sistema a ser modelado. Isso é muito vantajoso, pois, há uma menor possibilidade
de incluir agrupamentos que insiram dinâmicas espúrias ao sistema. Também aumenta-se a possibilidade de indicar os agrupamentos corretos que expliquem as não linearidades do sistema e
descartar agrupamentos que o explicam, diminuindo assim a quantidade de termos do modelo e
consequentemente simplif cando sua estrutura.
6.1 Sugestões para trabalhos futuros
101
6.1 Sugestões para trabalhos futuros
Algumas propostas para trabalhos futuros são citadas:
1. Utilizar a técnica apresentada neste trabalho como auxílio à outras técnicas de seleção de
estrutura;
2. Aplicar um sinal de frequências múltiplas na entrada para avaliação do espectro de frequência da saída e assim determinar uma estrutura de um modelo mais adequada;
3. Melhorar o cálculo de ErroW procurando uma melhor aproximação das frequências fundamental e ressonantes apresentadas pelo modelo em relação as frequências apresentadas
pelo sistema;
4. Inserir a análise de fase da resposta em frequência do sistema e dos agrupamentos candidatos relacionando-os com os atrasos máximos de entrada e saída.
102
Referências Bibliográf cas
Referências Bibliográficas
AGUIRRE, L. A. On the structure of nonlinear polynomial models: higher order correlation functions, spectra and term clusters. Circuits and Systems I - Fundamental Theory and
Applications, 1997. 450-453p. v.44.
AGUIRRE, L. A. Introdução à identificação de sistemas técnicas lineares e não-lineares
aplicadas a sistemas reais. UFMG, 2004.
AGUIRRE, L. A.; BILLINGS, S. A. Improved structure selection for nonlinear models based
on term clustering. International Journal of Control, v.62, n.3, p.569–587, 1995.
AGUIRRE, L. A.; BILLINGS, S. A. Dynamical effects of overparametrization in non-linear
models. Physica D, v.80, n.1-2, p.26–40, January 1995.
AGUIRRE, L. A.; BILLINGS, S. A. Improved structure selection for nonlinear models based
on term clustering. International Journal of Control, v.62, n.3, p.569–587, January 1995b.
AKAIKE, H. A new look at the statistical model identif cation. Automatic Control, IEEE
Transactions on, v.19, n.6, p.716–723, 1974.
AMARAL, G. F. V. Uso de redes neurais e conhecimento a priori na identificação de
sistemas dinâmicos não-lineares. 2001. — Programa de Pós-Graduação em Engenharia
Elétrica, Universidade Federal de Minas Gerais, Belo Horizonte, Brasil.
BILLINGS, S. A. Identif cation of nonlinear systems - a survey. IEE Proceedings, Parte D,
v.127, n.6, p.272–285, November 1980.
BILLINGS, S. A.; AGUIRRE, L. A. Effects of the sampling time and identif cation of nonlinear
models. International Journal of Bifurcation and Chaos, v.5, n.6, p.1541–1556, 1995.
BILLINGS, S. A.; CHEN, S. Representation of non-linear systems: the NARMAX model.
International Journal of Control, v.49, n.3, p.1013–1032, 1989.
104
Referências Bibliográf cas
BILLINGS, S. A.; CHEN, S. Identif cation of non-linear rational systems using a predictionerror estimation algorithm. International Journal of Systems Science, v.20, n.3, p.467–494,
1989.
BILLINGS, S. A.; CHEN, S.; KORENBERG, M. J. Identif cation of MIMO non-linear systems
using a forward-regression orthogonal estimator. International Journal of Control, v.49,
n.6, p.2157–2189, 1989.
BILLINGS, S. A.; TAO, Q. H. Model validity testes for non-linear signal processing applications. International Journal of Control, v.54, n.1, p.157–194, 1991.
BILLINGS, S. A.; ZHU, Q. M. Nonlinear Model Validation Using Correlations Tests.
Sheff eld S1 4DU, UK: University of Sheff eld, 1993. (463).
BILLINGS, S. A.; ZHU, Q. M. Nonlinear Model Validation using Correlation Tests. International Journal of Control, v.60, n.6, p.1107–1120, 1994.
BOYD, S.; CHUA, L. O.; DESOER, C. A. Analytical Foundations of Volterra Series. IMA
Journal of Mathematical and Information, v.1, p.243–282, 1984.
BRAGA, A. P.; CARVALHO, A. P. L. F.; LUDEMIR, T. B. Redes neurais artificiais: teoria e
aplicações. LTC, 2000.
CAMPOS, R. C. C. Projeto e construção de planta piloto de neutralização de pH e proposta
de metodologia para incorporação de informações auxiliares na identificação NARX
racional. 2007. — Programa de Pós-Graduação em Engenharia Industrial, Centro Universitário do Leste de Minas Gerais - Unileste-MG, Coronel Fabriciano, Brasil.
CHEN, S.; BILLINGS, S. A.; LUO, W. Orthogonal least squares methods and their applications
to non-linear system identif cation. International Journal of Control, v.50, n.5, p.1873–
1896, 1989.
COELHO, M. C. S. Modelos de Hammertein e de Wiener: conexões com modelos NARX
e sua aplicação em identif cação de sistemas não-lineares. 2002. — Programa de PósGraduação em Engenharia Elétrica, Universidade Federal de Minas Gerais, Belo Horizonte,
Brasil.
CORRÊA, M. V. Identificação caixa-cinza de sistemas não-lineares utilizando representações NARMAX racionais e polinomiais. 2001. — Programa de Pós-Graduação em
Engenharia Elétrica, Universidade Federal de Minas Gerais, Belo Horizonte, Brasil. Programa de Pós-Graduação em Engenharia Elétrica.
Referências Bibliográf cas
105
CORRÊA, M. V.; AGUIRRE, L. A. Identificação não-linear caixa cinza: uma revisão e novos
resultados. Revista Controle & Automação, 2004. (2, v.15).
GEVERS, M. Identif cation for control: from the early achievements to the revival of experiment design. European Journal of Control, v.11, n.4-5, p.335–352, 2005.
HABER, R.; UNBEHAUEN, H. Structure identif cation of nonlinear dynamic systems - a survey on input/output approaches. Automatica, v.26, n.4, p.651–677, 1990.
HAYKIN, S. Redes Neurais: princípios e prática. 2a .ed. Porto Alegre: Bookman, 2001.
HAYKIN, S.; VAN VEEN, B. Sinais e Sistemas. Porto Alegre: Bookman, 2001.
JING, X. J.; LANG, Z. Q.; BILLINGS, S. A. Output Frequency Response Function Based
Analysis for Nonlinear Volterra Systems. Mechanical Systems and Signal Processing, v.22,
n.1, p.102–120, 2008.
KORENBERG, M.; BILLINGS, S. A.; LIU, Y. P.; MCILROY, P. J. Orthogonal parameter estimation algorithm for non-linear stochastic systems. International Journal of Control, v.48,
n.1, p.193–210, 1988.
KORENBERG, M. J. Orthogonal identif cation of nonlinear difference equation models. Proceedings of the 28th Midwest Symposium on Circuits and Systems, p.90–95, 1985.
KORENBERG, M. J. A fast orthogonal search method for biological time-series analysis and
system identif cation. IEEE Inter. Conf. on Systems, Man and Cybernetics, v.1, n.1,
p.459–465, 1989a.
KORENBERG, M. J. A robust orthogonal algorithm for system identif cation and time-series
analysis. Biological Cybern. 1, v.69, n.1, p.267–276, 1989b.
LANG, Z. Q.; BILLINGS, S. A. Evaluation of Output Frequency Responses of Nonlinear Systems under Multiple Inputs. IEEE Trans. on Circuits and Systems, v.47, n.1, p.28–38,
2000.
LANG, Z. Q.; BILLINGS, S. A.; YUE, R.; LI, J. Output Frequency Response Function of
Nonlinear Volterra Systems. Automatica, v.43, p.805–816, 2007.
LEONTARITIS, I. J.; BILLINGS, S. A. Input-output parametric models for nonlinear systems Part I: deterministic non-linear systems. International Journal of Control, v.41, n.2, p.303–
328, 1985.
106
Referências Bibliográf cas
LEONTARITIS, I. J.; BILLINGS, S. A. Input-output parametric models for nonlinear systems
- Part II: stochastic nonlinear systems. International Journal of Control, v.41, n.2, p.329–
344, 1985.
LEONTARITS, I. J.; BILLINGS, S. A. Model selection and validation methods for non-linear
systems. International Journal of Control, v.45, n.1, p.311–341, 1987.
LJUNG, L. System Identification: theory for the user. 2.ed. Englewood Cliffs, New Jersey:
Prentice-Hall, 1999.
LU, S.; JU, K. H.; CHON, K. H. A new algorithm for linear and nonlinear ARMA model
parameter estimation using aff ne geometry. IEEE Trans. on Biomedical Engineering, v.48,
n.10, p.1116–1124, 2001.
MAO, K. Z.; BILLINGS, S. A. Algorithms for minimal model structure detection in nonlinear
dynamic system identif cation. International Journal of Control, v.68, p.311–330, 1997.
MATKO, D.; ZUPANCIC, B.; KARBA, R. Simulation and modeling of continuous systems
- A case study approach. Prentice Hall International, 1992.
MENDES, E. M. A. M.; BILLINGS, S. A. An alternative solution to the model structure selection problem. IEEE Transaction on Systems, Man, and Cybernetics - Part A: Systems
and Humans, v.31, n.6, p.597–608, November 2001.
OGATA, K.; KOHN, A. F.; BARROS MORAIS, J. C. T. de. Engenharia de controle moderno.
Prentice Hall, 2003.
PAARMANN, L. D.; KORENBERG, M. J. Estimation of the parameters of an ARMA signal
model based on an orthogonal search. IEEE Trans. on Automatic Control, v.37, n.3, p.347–
352, 1992.
PEARSON, R.; POTTMANN, M. Gray-box identif cation of block oriented nonlinear models.
Journal of Process Control, v.10, n.4, p.301–315, August 2000.
PEYTON-JONES, J. E.; BILLINGS., S. Recursive algorithm for computing the frequency response of a class of non-linear difference equations models. International Journal of Control, v.50, n.5, p.1925–1940, 1989.
PIRODDI, L. Simulation error minimisation methods for NARX model identif cation. International Journal of Modeling, Identification and Control, v.3, n.4, p.392–403, 2008.
Referências Bibliográf cas
107
PIRODDI, L.; SPINELLI, W. An identif cation algorithm for polynomial NARX models based
on simulation error minimization. International Journal of Control, v.76, n.17, p.1767–
1781, 2003.
RHODES, C.; MORARI, M. Determining the model order of nonlinear input/output systems
directly from data. Proceedings of the American Control Conference, v.3, p.2190–2194,
1995.
RHODES, C.; MORARI, M. Determination of model order for NARX models directly from
input-output data. Journal of Process Control, v.8, n.5,6, p.459–468, 1995.
RODRÍGUEZ-VÁZQUEZ, K.; FONSECA, C. M.; ; FLEMING, P. J. Identifying the structure of nonlinear dynamic systems using multiobjective genetic programming. IEEE Transactions on Systems, Man and Cybernetics - Part A: Systems and Humans, v.34, n.4,
p.531–545, 2004.
RUGH, W. J. Nonlinear systems theory - the Volterra/Wiener approach. 1981.
SJOBERG, J.; ZHANG, Q.; LJUNG, L.; BENVENISTE, A.; DELYON, B.; GLORENNEC,
P.; HJALMARSSON, H.; JUDITSKY, A. Nonlinear black-box modeling in system identif cation: a unif ed overview. Automatica, v.31, n.12, p.1691–1724, 1995.
SPINELLI W., P.-L. e. L. M. A two-stage algorithm for structure identif cation of polynomial
NARX models. In: 2006, Minneapolis - Minnesota. Anais. . . 2006. p.2387–2392.
TULLEKEN, H. J. A. F. Grey-box modelling and identif cation using physical knowledge and
Bayesian techniques. Automatica, Tarrytown, NY, USA, v.29, n.2, p.285–308, 1993.
VOLTERRA, V. Theory of functions. Blackie and Sons, 1930.
WEI, H. L.; BILLINGS, S. A. Model structure selection using an integrated forward orthogonal search algorithm assisted by squared correlation and mutual information. International
Journal of Modeling, Identification and Control, v.4, n.3, p.341–346, 2008.
WEI H. L. E BILLINGS, S. A. Model structure selection using an integrated forward orthogonal search algorithm assisted by squared correlation and mutual information. International
Journal of Modeling, Identification and Control, v.3, n.4, p.341–356, 2008.
M.I.T., T. T. P. (Ed.). Nonlinear problems in random theory. New York, NY.: In Wiley, J. and
Sons, editors, 1958.
108
Referências Bibliográf cas
WIGREN, T. Recursive prediction error identif cation using nonlinear Wiener model. Automatica, v.29, n.4, p.1011–1025, 1993.
WU, X.; LANG, Z. Q.; BILLINGS, S. A. Analysis of Output Frequencies of Nonlinear Systems. IEEE Trans. on Signal Processing, v.55, n.7, p.3239–3246, 2007.
ZHANG L. F., Z. Q. M. e. L. A. A set of novel correlation tests for nonlinear system variables.
International Journal of Systems Science, v.38, n.1, p.47–60, 2007.
ZHANG, L. F.; ZHU, Q. M.; LONGDEN, A. A set of novel correlation tests for nonlinear
system variables. International Journal of Systems Science, v.38, n.1, p.47–60, 2007.
ZOU, R.; CHON, K. H. Robust algorithm for estimation of time-varying transfer functions.
IEEE Trans. on Biomedical Engineering, v.51, n.2, p.219–228, 2004.