ISSN 1984-8218
APLICAÇÃO DO ALGORITMO DOS VAGALUMES NA
IDENTIFICAÇÃO SIMULTÂNEA DA ESPESSURA ÓPTICA E
ALBEDO COM VARIAÇÃO ESPACIAL EM UM PROBLEMA
INVERSO DE TRANSFERÊNCIA RADIATIVA
Rubens L. Cirino
Instituto Politécnico – IPRJ – Nova Friburgo
Universidade do Estado do Rio de Janeiro – UERJ
Diego C. Knupp
Agência Nacional de Transportes Terrestres – ANTT
Antônio J. Silva Neto
Instituto Politécnico – IPRJ – Nova Friburgo
Universidade do Estado do Rio de Janeiro - UERJ
Resumo: A necessidade de solução de problemas de otimização aparece em diferentes
contextos com aplicações práticas, notadamente em engenharia, onde tem importância no
contexto de otimização de desenho de produtos, layout de produção e logística. Outro contexto
onde a solução de problemas de otimização exerce papel fundamental é na solução de
problemas inversos, quando formulados implicitamente e a solução é dada pela minimização de
uma função objetivo dada pelo somatório dos resíduos entre as medidas experimentais e as
grandezas calculadas em função dos parâmetros sendo estimados. Neste trabalho, motivados
por um problema inverso em transferência radiativa, introduzimos a utilização da heurística
conhecida como algoritmo dos vagalumes (Firefly Algorithm – FA) na minimização da função
objetivo para solução do problema inverso de estimativa simultânea da espessura óptica e do
albedo com variação espacial, formulado como um problema de estimativa de parâmetros. Os
resultados são criticamente comparados a dois outros métodos estocásticos de otimização
bastante difundidos na comunidade científica, o algoritmo de colisões de partículas (Particle
Collision Algorithm – PCA) e o algoritmo de enxame de partículas (Particle Swarm
Optimization – PSO), que já foram utilizados com sucesso pelos autores na solução de
problemas inversos em transferência radiativa. O desempenho superior do algoritmo dos
vagalumes observado no caso teste apresentado motiva o estudo mais aprofundado,
hibridizações e derivações com base nesta heurística em trabalhos futuros.
Introdução
Os problemas inversos são geralmente formulados implicitamente [1], definindo-se uma função
objetivo dada pela soma dos quadrados dos resíduos entre as medidas experimentais e os
correspondentes valores calculados através do modelo direto, em função dos parâmetros que se
deseja estimar. Assim, a solução do problema inverso se torna o problema de minimização da
função objetivo definida. Dada a forma funcional da função objetivo, geralmente complexa e
repleta de mínimos locais, o emprego direto de métodos determinísticos de otimização pode
levar à não convergência da solução ou à estagnação em um dos mínimos locais. Métodos
estocásticos de otimização, geralmente inspirados em comportamentos observados na natureza,
têm se apresentado eficazes na recuperação de mínimos globais [2], mesmo de funções
complexas, dado um número suficientemente alto de iterações. O desenvolvimento paulatino na
velocidade de processamento dos computadores vem permitindo tanto a utilização de métodos
estocásticos de otimização quanto a proposição de modelos computacionais cada vez mais
complexos em problemas de engenharia, motivando o desenvolvimento de novas heurísticas e
suas aplicações em uma vasta gama de problemas.
O estudo de problemas inversos aplicados à transferência radiativa tem sido objeto de intensa
pesquisa objetivando atender a diversas aplicações que vão desde a medicina até a indústria [3],
1309
ISSN 1984-8218
onde geralmente se busca a realização de testes não invasivos/destrutivos [4]. Nosso grupo de
pesquisa publicou diversos trabalhos em transferência radiativa onde buscou-se comparar
métodos estocásticos de otimização já estabelecidos, bem como testar novas propostas de
hibridização e variações de métodos tradicionais. Apenas para citar alguns exemplos,
abordamos o PCA e o Luus-Jaakola [5], o PSO [6], O recozimento simulado (Simulated
Annealing – SA), algoritmos genéticos (Genetic Algorithms – GA), otimização por colônia de
formigas (Ant Colony Optimization – ACO), otimização extrema generalizada (Generalized
Extremal Optimization – GEO) e evolução diferencial (Differential Evolution – DE) [2].
O presente trabalho está focado no estudo do Firefly Algorithm (FA) para a solução de um
problema em transferência radiativa onde é considerado que o albedo de espalhamento
apresenta variação espacial e deve ser estimado [7]. Neste trabalho também consideramos a
estimativa simultânea da espessura óptica do meio [8]. Com o intuito de ter uma medida
objetiva do desempenho do FA, comparamos suas soluções com dois outros métodos
estocásticos já bastante difundidos na comunidade científica, PCA e PSO, cuja aplicação no
tratamento de problemas inversos em transferência radiativa já foi abordada com sucesso [2,6].
Formulação Matemática e Solução do Problema Direto
Considera-se um meio unidimensional, cinza, heterogêneo, com espalhamento isotrópico de
espessura óptica τ0 , com superfícies de contorno transparentes. Os contornos em τ = 0 e τ = τ0
estão sujeitos a fontes externas de radiação isotrópica com intensidades A1 e A2,
respectivamente, como ilustra esquematicamente a Fig. 1, onde ρ1 = ρ 2 = 0 .
Figura 1. Representação esquemática de um meio participante unidimensional sujeito à
incidência de radiação originada por fontes externas.
O modelo matemático para a interação da radiação com o meio participante é dado pela equação
de Boltzmann, que para o caso de simetria azimutal e albedo com dependência espacial é escrito
na forma adimensional como:
∂I (τ , u )
ω ( x) 1
+ I (τ , u ) =
I (τ , u ' )dµ ' ;
∫
−
1
∂τ
2
I (0, µ ) = A1 , µ > 0, I (τ 0 , µ ) = A2 , µ < 0
µ
em 0 < τ < τ0 , para -1≤µ≤ 1
(1a)
(1b,c)
onde τ é a variável óptica, I é a intensidade da radiação, µ é o cosseno do ângulo polar e ω (τ ) é
o albedo de espalhamento que aqui é representado como a seguinte expansão:
K
ω (τ ) = ∑ Dkτ k
(2)
k =0
Quando a geometria, as propriedades radiativas e as condições de contorno são conhecidas, o
problema (1) pode ser resolvido e a intensidade da radiação I é determinada para todo o domínio
espacial e angular, isto é, 0 ≤ τ ≤ τ0 e -1 ≤ µ ≤ 1. Este é o chamado problema direto. Quando as
1310
ISSN 1984-8218
propriedades radiativas ou as condições de contorno não são conhecidas, mas dados
experimentais das intensidades de radiação podem ser obtidos, os parâmetros desconhecidos
podem ser estimados através da proposição de um problema inverso.
Formulação Matemática e Solução do Problema Inverso
Considere a espessura óptica do meio, τ0 , e o albedo com dependência espacial,
representado na forma dada pela Eq. (2), sejam desconhecidos. Ou seja, desejamos encontrar
, k = 0,1,2,..., K. Portanto, temos o seguinte vetor de
estimativas para τ0 e os coeficientes
incógnitas:
Considere que os dados experimentais adquiridos em ambos os contornos, externamente ao
meio, em diferentes ângulos polares, podem ser obtidos. Assim, temos a disposição Yi, i = 1, 2,
3,..., Nd, onde Nd é o número total de dados experimentais. Este problema inverso pode, então,
ser formulado como um problema de otimização, onde buscamos minimizar a função objetivo
dada pelo somatório dos quadrados dos resíduos:
Neste trabalho, dados experimentais reais não estão disponíveis e as intensidades experimentais,
, são então simuladas, calculando-se os valores com a solução do problema (1) usando os
valores exatos dos parâmetros que desejamos estimar no problema inverso. A estes valores são
adicionados ruídos aleatórios de uma distribuição normal com média zero, simulando erros
experimentais:
(5)
onde r são números aleatórios de uma distribuição normal com média zero e desvio padrão
unitário, e σ e simula o desvio padrão dos erros de medição.
Para a minimização da função objetivo dada pela Eq. (4), neste trabalho utilizamos o algoritmo
dos vagalumes, cuja descrição é dada em maiores detalhes a seguir. De modo a ser objetivo
quanto ao desempenho deste método na solução do problema investigado, fazemos uma
comparação crítica com dois outros métodos estocásticos já bastante difundidos, o algoritmo de
colisões de partículas (Particle Collision Algorithm – PCA) e o algoritmo de enxame de
partículas (Particle Swarm Optimization – PSO). Por questões de brevidade, as descrições
destes dois últimos métodos serão omitidas, mas podem ser encontradas em detalhes nas
referências [9-12].
Algoritmo dos Vagalumes (Firefly Algorithm - FA)
Para uma melhor compreensão do Firefly Algorithm (FA) [11,12], duas características do
algoritmo devem ser destacadas: como se dá a variação da intensidade da luz percebida pelo
vagalume; e como é formulada a atratividade entre os vagalumes. A intensidade de emissão de
luz por parte de um vagalume é proporcional à função objetivo, porém a intensidade de luz
percebida por um vagalume decai em função da distância entre os vagalumes. Logo, a
intensidade percebida por um vagalume é dada por:
, em que é a intensidade
da luz emitida; r é a distância Euclidiana entre os vagalumes i e j, sendo i o vagalume mais
brilhante e j o vagalume menos brilhante; γ é o parâmetro de absorção da luz pelo meio.
Desta maneira o fator de atratividade pode ser formulado como:
(6)
1311
ISSN 1984-8218
onde
é a atratividade para uma distância r = 0, e pode ser fixado em
. Assim, a
movimentação em um dado passo de tempo t de um vagalume i em direção a um melhor
vagalume j é definida como:
onde o segundo termo do lado direito da equação insere o fator de atratividade enquanto o
terceiro termo, regulado pelo parâmetro , insere aleatoriedade no caminho percorrido pelo
vagalume; rand é um número aleatório de uma distribuição uniforme entre 0 e 1. O
Pseudocódigo do Firefly Algorithm é apresentado na Fig. 2.
Início
Definir a função objetivo;
Definir os parâmetros
;
Gerar a população inicial de vagalumes
Para t=1 até
;
Calcular a intensidade para
proporcionalmente à função
objetivo;
Para i=1 até n;
Calcular o fator de atratividade de acordo com
;
Mover o vagalume i em direção aos vagalumes mais
brilhantes.
Verificar se o vagalume está dentro dos limites;
Fim-Para
Fim-Para
Pós-processar e visualizar os resultados
Fim
Figura 2. Pseudocódigo do Firefly Algorithm.
Resultados e Discussão
Nos resultados apresentados para o caso teste selecionado a seguir, foram utilizados dois níveis
e
,
de ruído diferentes na simulação dos dados experimentais, com
resultando em erros de até 2% e 5%, respectivamente. O problema selecionado para os testes
possui espessura óptica unitária, isto é, τ 0 = 1 e o albedo com variação espacial é dado por uma
reta, ω (τ ) = τ . Consideramos ainda a iluminação externa dada por A1 = 1 e A2 = 0 , nas Eqs.
(1b,c).
No procedimento de solução do problema inverso, consideramos que a expansão dada pela Eq.
(2) possui três termos, ou seja, K = 3. Portanto, para o exemplo investigado, temos os seguintes
valores exatos dos parâmetros a serem estimados:
r
Z exato = {τ 0 = 1, D0 = 0, D1 = 1, D2 = 0}
(8)
Na solução do problema de otimização, consideramos os seguintes parâmetros para o FA:
e γ = 1 . A escolha destes parâmetros foi baseada em
n = 20 , MaxGerações = 20 ,
[11]. O PSO, por sua vez, foi ajustado considerando uma população de 80 partículas por 20
gerações e para os parâmetros foi considerado α = β = 2 e γ = 1 , escolhidos também com
base na ref. [11]. Quanto ao PCA, este método não apresenta a necessidade de ajuste de outros
parâmetros além do número de iterações no loop externo e interno, que aqui foram escolhidos
como 20 e 20, em ambos os casos. Os parâmetros que definem o número de iterações foram
escolhidos de modo que os três métodos desempenhassem aproximadamente o mesmo esforço
computacional, neste caso, em torno de 4500 avaliações da função objetivo. Os algoritmos
foram executados em um PC com processador Intel Core2 Duo CPU T6670 @2.20GHz,
2.96GB RAM, e sistema operacional Windows 7.0 de 32 bits. Cada execução, para cada um dos
métodos, levou aproximadamente 23 minutos.
1312
ISSN 1984-8218
Para cada método, a solução do problema de otimização foi executada por 11 vezes de forma
independente, considerando o mesmo conjunto de dados experimentais. Como os métodos são
estocásticos, esperamos uma solução distinta a cada execução e o objetivo de várias execuções é
medir a dispersão das soluções obtidas. Nas Tabs. 1 e 2 a seguir, apresentamos o resumo das
soluções para os três métodos estudados, considerando 2 e 5% de ruído nos dados
e
se referem à melhor e à pior
experimentais, respectivamente. Nestas tabelas,
estimativa obtida dentre as 11 execuções, respectivamente. Aqui, melhores estimativas são
consideradas quando levam a um menor valor na avaliação da função objetivo.
Nestes resultados, fica evidente que o FA foi o que apresentou melhor desempenho na
minimização da função objetivo entre os três algoritmos investigados, para ambos os níveis de
ruído, onde suas melhores estimativas levaram a função objetivo a valores até três ordens de
grandeza inferiores ao PCA e o PSO. Outra importante conclusão que pode ser inferida destas
tabelas diz respeito à dispersão das soluções. O desvio padrão das estimativas obtidas para a
espessura óptica, τ 0 , onde os três métodos tiveram bom desempenho – resultando em
estimativas bastante próximas do esperado – o FA demonstrou menor dispersão nas soluções
dentre os três casos. Quanto à dispersão nas estimativas dos coeficientes D0 , D1 e D2 , os
desvios padrão calculados são maiores para o FA com relação ao PCA, entretanto, tanto a média
quanto a melhor solução do FA são muito mais próximas do valor esperado, indicando que o
PCA provavelmente ficou estagnado na redondeza do mesmo mínimo local em várias das 11
execuções.
Ressalta-se ainda, que quanto aos coeficientes D0 , D1 e D2 , o objetivo final é a estimativa da
função ω (τ ) a partir destes coeficientes, e não a estimativa dos parâmetros em si. Portanto, as
Figs. 3a,b apresentam as curvas traçadas a partir das estimadas na melhor execução de cada
método em comparação com a curva exata, ω (τ ) = τ , para dados experimentais simulados com
2% e 5% de erro, respectivamente, onde fica evidente que os três métodos foram capazes de
obter estimativas para os coeficientes que geram curvas que tendem para a curva exata.
Entretanto fica claro que a curva estimada pelo FA foi a que ficou mais próxima do resultado
esperado, praticamente coincidente com a curva exata.
Finalmente, as Figs. 4a,b apresentam a evolução do valor da função objetivo, dada pela solução
naquele momento, com respeito ao número de acessos. Estes resultados confirmam que o FA
apresentou desempenho superior ao PCA e ao PSO, que parecem ficar estagnados em uma
solução por um grande número de iterações, provavelmente um mínimo local. Cabe ressaltar
que o número reduzido de iterações teve justamente o objetivo de evidenciar diferenças de
desempenho entre os três métodos. Outro ponto a ser ressaltado, é o melhor desempenho do FA
alcançado com um ruído de 5% e não com o de 2%. Este resultado não era o esperado, o que
merece, ainda, um estudo mais aprofundado para se tentar levantar a sua causa.
Conclusões
Neste trabalho, motivados por um problema inverso em transferência radiativa, introduzimos a
utilização da heurística conhecida como algoritmo dos vagalumes (Firefly Algorithm – FA) na
minimização da função objetivo para solução do problema inverso de estimativa simultânea da
espessura óptica e do albedo de espalhamento com variação espacial, formulado como um
problema de estimativa de parâmetros. Os resultados são criticamente comparados a dois outros
métodos estocásticos de otimização bastante difundidos na comunidade científica, o algoritmo
de colisões de partículas (Particle Collision Algorithm – PCA) e o algoritmo de enxame de
partículas (Particle Swarm Optimization – PSO), que já foram utilizados com sucesso pelos
autores na solução de problemas inversos em transferência radiativa. Os resultados obtidos pelo
FA são notoriamente superiores aos obtidos pelo PCA e PSO, indicando boa perspectiva na
utilização desta heurística nesta classe de problemas. O trabalho deve prosseguir na investigação
1313
ISSN 1984-8218
mais aprofundada dos mecanismos de inteligência artificial do FA para desenvolvimento de
variações e hibridizações específicas, de modo a garantir sua robustez e desempenho,
características desejáveis na solução de problemas computacionalmente intensivos.
(a)
(b)
Figura 3. Estimativas de ω (τ ) em comparação com a curva exata para ruídos de: (a) até 2%.
(b) até 5%.
(a)
(b)
Tabela 1. Resultados obtidos com o FA, PCA e PSO para dados experimentais simulados com
ruído de: (a) até 2% (b) até 5%.
Agradecimentos
Os autores agradecem ao suporte financeiro do Conselho Nacional de Desenvolvimento
Científico e Tecnológico, CNPq e da Fundação Carlos Chagas Filho de Amparo à Pesquisa do
Estado do Rio de Janeiro, FAPERJ.
1314
ISSN 1984-8218
(a)
(b)
Figura 4. Evolução do valor da função objetivo para níveis de ruído de: (a) até 2%. (b) até 5%.
Referências
[1] Silva Neto, A. J., e Moura Neto, F. D., Problemas Inversos – Conceitos Fundamentos e
Aplicações, Ed. UERJ, 2005.
[2] Silva Neto, A. J., e Becceneri, J. C., Técnicas de Inteligência Computacional Inspiradas na
Natureza – Aplicação em Problemas Inversos em Transferência Radiativa – Notas em
Matemática Aplicada – SBMAC – Vol. 41, 2009.
[3] Arridge, S.R., Optical tomography in medical imaging, inverse problems, Vol. 15, pp R41
– R93, 1999.
[4] Oliva Soares, P., Soeiro, F. J. C. P. e Silva Neto, A.J., Solution of inverse radiative transfer
problems with artificial neural networks and hybrid methods, Inverse Problems in
Engineering, Seminar, Cincinnati, USA, 2004.
[5] Knupp, D.C., Silva Neto, A.J., Sacco, W. F. – Radiative properties estimation with the
luus-jaakola and the particle collision algorithm, CMES, vol.54, no.2, pp.121-145, 2009.
[6] Becceneri, J. C., Stephany, S., Campos Velho, H. F., Silva Neto, A.J., Solution of the
inverse problem of radiative properties estimation with the particle swarm optimization
technique, 14th Inverse Problems In Engineering Seminar, 2006.
[7] Stephany, S, Becceneri, J. C., Souto, R. P., Campos Velho, H. F., Silva Neto, A. J., A preregularization scheme for the reconstruction of a spatial dependent scattering albedo using
a hybrid ant colony optimization implementation - Applied Mathematical Modelling 34pp. 561-572, 2010.
[8] Knupp, D.C., Silva Neto, A. J., Simultaneous identification of the optical thickness and
space-dependent albedo using Bayesian Inference. 21st International Congress of
Mechanical Engineering, 2011.
[9] Sacco, W. F. and de Oliveira, C. R. E. A new stochastic optimization algorithm based on a
particle collision metaheuristic. Proceedins of 6th WCSMO, 2005.
[10] Kennedy, J. and Eberhart, R. C. Particle swarm optimization. Proceedings of IEEE
International Conference on Neural Networks, Piscataway, NJ. pp. 1942-1948, 1995.
[11] Yang, Xin-She, Nature-Inspired Metaheuristic Algorithms, Luniver Press, 2008.
[12] Luz, E.F.P., Becceneri, J. C., Campos Velho, H. F., Conceitualização do algoritmo
vagalume e sua aplicação na estimativa de condição inicial da equação de calor; Workshop
dos Cursos de Computação Aplicada, São José dos Campos/SP – Brasil, 2009.
1315