FABIANNO GIANELLI DELALIBERA Aplicação do algoritmo Artificial Bee Colony na otimização inversa de cristais fotônicos bidimensionais São João da Boa vista - SP 2020 FABIANNO GIANELLI DELALIBERA Aplicação do algoritmo Artificial Bee Colony na otimização inversa de cristais fotônicos bidimensionais Dissertação apresentada como parte dos requisitos para obtenção do título de Mestre em Engenharia Elétrica, junto ao Programa de Pós- Graduação em Engenharia Elétrica, Interunidades, entre o Instituto de Ciência e Tecnologia de Sorocaba e o Câmpus de São João da Boa Vista da Universidade Estadual Paulista “Júlio de Mesquita Filho”. Orientador: Prof. Dr. Gilliard Nardel Malheiros Silveira São João da Boa Vista - SP 2020 Sistema de geração automática de fichas catalográficas da Unesp. Biblioteca do Câmpus Experimental de São João da Boa Vista. Dados fornecidos pelo autor(a). Essa ficha não pode ser modificada. D335a Delalibera, Fabianno Gianelli Aplicação do algoritmo Artificial Bee Colony na otimização inversa de cristais fotônicos bidimensionais / Fabianno Gianelli Delalibera. -- São João da Boa Vista, 2020 71 p. : il., tabs. Dissertação (mestrado) - Universidade Estadual Paulista (Unesp), Câmpus Experimental de São João da Boa Vista, São João da Boa Vista Orientador: Gilliard Nardel Malheiros-Silveira 1. Otimização matemática. 2. Fotônica. 3. Ondas eletromagneticas. I. Título. FABIANNO GIANELLI DELALIBERA Aplicação do algoritmo Artificial Bee Colony na otimização inversa de cristais fotônicos bidimensionais Dissertação apresentada como parte dos requisitos para obtenção do título de Mestre em Engenharia Elétrica, junto ao Programa de Pós- Graduação em Engenharia Elétrica, Interunidades, entre o Instituto de Ciência e Tecnologia de Sorocaba e o Câmpus de São João da Boa Vista da Universidade Estadual Paulista “Júlio de Mesquita Filho”. Comissão Examinadora Prof. Dr. Gilliard Nardel Malheiros Silveira UNESP – Câmpus de São João da Boa Vista. Orientador Prof. Dr. Daniel Orquiza de Carvalho UnB – Universidade de Brasília. Prof. Dr. Vitaly Félix Rodríguez Esquerre UFBA – Universidade Federal da Bahia. São João da Boa Vista - SP 23 de julho de 2020 “E a perseverança deve ter ação completa, a fim de que vocês sejam maduros e íntegros, sem que falte a vocês coisa alguma.” Tiago 1:4 Bíblia Sagrada AGRADECIMENTOS Agradeço primeiramente a Deus, pelas bençãos concedidas a mim no decorrer de toda a minha vida. Em especial, por ter providenciado que o Prof. Dr. Gilliard Nardel Malheiros Silveira tenha me aceitado como seu orientando. Trabalhar ao seu lado me permitiu admirar e respeitar o seu trabalho e dedicação, professor. Portanto, ao Prof. Dr. Gilliard Nardel Malheiros Silveira, minha gratidão pela oportuni- dade de ser seu orientando, por permitir conciliar estudo e trabalho, garantindo uma excelente gestão de tempo, e também pelos conselhos oferecidos ao longo do curso. Meu agradecimento sincero ao Prof. Dr. Vitaly F. Rodríguez-Esquere por disponibilizar o software de cálculo eletromagnético baseado no FEM. A minha esposa, Suellen, agradeço pela paciência e por se privar da minha presença aos fins de semana, me permitindo dedicar aos estudos. Obrigada por sempre me apoiar e incentivar. Aos meus pais, José Roberto e Genoefa, e minha irmã Daniella, muito obrigada por todo apoio, carinho e amor oferecidos incondicionalmente. Vocês são uma inspiração para mim e meus exemplos de determinação. À UNIVESP, agradeço por oportunizar a experiência de atuar como facilitador em di- versas disciplinas no decorrer de um ano e meio, contando com o bem-vindo auxílio de uma bolsa de estudos. À todos professores com os quais tive o privilégio de trabalhar, absorvendo um pouco de seus conhecimentos generosamente compartilhados em cada disciplina cursada: meu muito obrigado! Agradeço também ao meu professor de inglês, Geoclândio B. dos Santos Filho, pelo auxílio e ensinamentos da língua Inglesa. E por fim, agradeço à toda a diretoria da Dacota Condutores Elétricos - Ltda e, em especial, ao meu gerente Eng. Luis M. Machado, por disponibilizarem a flexibilidade de horário de trabalho, me permitindo participar das disciplinas do curso mestrado. 8 RESUMO Neste trabalho, apresenta-se uma nova proposta de aplicação do algoritmo Artificial Bee Colony (ABC) em um problema eletromagnético para otimização de estruturas cristalinas fotônicas bidimensionais em arranjo triangular, buscando maximizar a banda proibida fotônica em estruturas compostas por silício e ar, tanto para o modo transversal elétrico quanto para o transversal magnético. Estruturas otimizadas via Genetic Algorithm (GA), apresentadas na literatura, foram utilizadas como parâmetro comparativo. O estudo foi realizado, inicialmente, submetendo os algoritmos clássicos do ABC e GA à busca pelo valor máximo de duas funções bidimensionais previamente conhecidas, com o objetivo de compreender o comportamento de busca do ABC para, posteriormente, ser configurado à otimização de estruturas fotônicas. Foram necessárias algumas modificações no ABC clássico, para otimização de problemas em um espaço de buscas discreto. Os resultados apresentados pelo ABC foram bastante satisfatórios, tendo encontrado estruturas com percentual de banda proibida fotônica superiores ao encontrado pelo GA. Sendo assim, o ABC apresentou uma melhor eficiência para o problema proposto em comparação ao GA. Palavras-chaves: Algoritmo Colônia de Abelhas Artificiais, Cristal Fotônico, Banda Proibida Fotônica. 9 ABSTRACT In this work, it is show a novel proposal for the application of the Artificial Bee Colony (ABC) algorithm in an electromagnetic problem, for optimization of two-dimensional photonic crystal structure in triangular lattice, seeking to maximize the photonic band gap (PBG) in structures composed of silicon and air, both for the modes transverse electric and transverse magnetic. Structures optimized by Genetic Algorithm (GA), presented in the literature, were used as a comparative parameter. The study was initially carried out by submitting the classical algorithms of ABC and GA to search for the maximum value of two previously known two-dimensional functions, with the objective of understanding the search behavior of ABC to later configure the optimization of photonic structures. Some modifications were needed in the classic ABC to optimize problems in a discrete search space. The results presented by the ABC were quite satisfactory, finding structures with a percentage of PBG higher than that found by GA. Therefore, the ABC presented a better efficiency for the proposed problem compared to GA. Keywords: Artificial Bee Colony, Photonic Crystal, Photonic Band Gap. 10 LISTA DE ILUSTRAÇÕES Figura 1 – Exemplo de operação crossover entre dois indivíduos, com apenas um ponto de crossover e um bit de mutação. .......................................................... 23 Figura 2 – Porcentagem de publicação relacionada a cada Algoritmo HBMO, BA, BCO e ABC. ....................................................................................................... 26 Figura 3 – Publicações em diferentes áreas nos cinco primeiros anos do algoritmo ABC. ................................................................................................................... 27 Figura 4 – Fluxograma do comportamento do ABC separado em etapas. ................ 29 Figura 5 – Exemplo de cristais fotônicos unidimensional, bidimensional e tridimensional. .................................................................................................... 33 Figura 6 – (a) Estrutura cristalina com periodicidade discreta Λ na direção �̂� e (b) propagação do campo 𝐸(x,z) entre a estruturas para os três principais regimes: sub-comprimento de ondas, reflexão e difração. ............................................... 34 Figura 7 – Estrutura de uma cristal fotônico em duas dimensões com arranjo quadrado (a) Rede de Bravais com seus vetores primitivos; (b) Os correspondentes pontos da rede recíproca e seus vetores primitivos; (c) Construção da primeira zona de Brillouin. ......................................................... 40 Figura 8 – Estrutura de banda fotônica de cristal fotônico unidimensional com a = 0,5Λ com relação entre as constante dielétrica de: (a) Δɛ = 1; (b) Δɛ = 13/12; (c) Δɛ = 13/1. ........................................................................................................... 42 Figura 9 – (a) Estrutura cristalina com arranjo quadrado de materiais dielétricos em formato de cilindro com raio r e periodicidade Λ, (b) Estrutura de banda fotônica para um arranjo quadrado com colunas de Alumina com n ≈ 2,98, envolto de ar. ........................................................................................................................... 43 Figura 10 - Guia de onda reto em uma estrutura cristalina fotônica bidimensional. (a) diagrama de banda para os vetores de onda na mesma direção da linha com defeito e (b) a propagação da luz com a polarização do campo �⃗� 𝑧. .................. 44 Figura 11 – Imagem microscopia de um guia de onda com curva de 90º em uma estrutura cristalina fotônica. ............................................................................... 44 Figura 12 – Exemplo de filtro passa faixa, (a) estrutura cristalina e suas respectivas relações geométricas (b) gráfico de transmissão em relação as frequências para relações: r1 = 0,29Λ e r2 = 0,18Λ. ....................................................................... 45 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107232 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107232 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107233 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107233 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107234 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107234 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107236 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107236 11 Figura 13 - Divisor de potência em cristal fotônico bidimensional. (a) Esquemático da estrutura formado por linhas de defeitos reduzindo o diâmetro das colunas de dielétrico; (b) Propagação do campo �⃗� 𝑧 dentre a estrutura. .............................. 46 Figura 14 – Gráfico de superfície Função 1. ............................................................. 48 Figura 15 – Resultados das simulações referente a Função 1 com todos os parâmetros de parada propostos para (a) 10 indivíduos, (b) 20 indivíduos, (c) 40 indivíduos e (d) 80 indivíduos. ............................................................................ 49 Figura 16 – Gráfico tridimensional da função 1 com as posições finais de cada indivíduo para 80 indivíduos e 800 ciclos. .......................................................... 50 Figura 17 – Gráfico de evolução para Função 1 com configuração de 80 indivíduos e 800 ciclos dos algoritmos. (a) Algoritmo ABC e (b) Algoritmo GA. ..................... 50 Figura 18 – Gráfico de superfície Função 2. ............................................................. 51 Figura 19 – Resultados das simulações para função 2 para cada configuração, dividindo em quatro gráficos para diferentes quantidades de população (a) 10 indivíduos, (b) 20 indivíduos, (c) 40 indivíduos e (d) 80 indivíduos. ................... 51 Figura 20 – Gráfico tridimensional da função 2 com as posições finais de cada indivíduo para 80 indivíduos e 800 ciclos. .......................................................... 52 Figura 21 – Gráfico de evolução para Função 1 com configuração de 80 indivíduos e 800 ciclos dos algoritmos. (a) Algoritmo ABC e (b) Algoritmo GA. ..................... 53 Figura 22 – Média do tempo computacional total utilizado pelo GA e o ABC, para configuração de 80 indivíduos e 800 ciclos, para: (a) Função 1 e (b) Função 2. ........................................................................................................................... 53 Figura 23 – Visão bidimensional de estruturas cristalinas fotônicas sendo formada por: (a) Lacunas cilíndricas de ar em arranjo quadrado; (b) Cilindros dielétricos em arranjo quadrado; (c) Lacunas cilíndricas de ar em arranjo triangular; (d) Cilindros dielétricos em arranjo triangular. ......................................................... 54 Figura 24 – (a) Estrutura com arranjo triangular cristalina composta por 16 células unitárias; (b) Célula unitária representado por 9 elementos e (c) Representação de uma abelha em modo binário para célula unitária ao lado. ........................... 55 Figura 25 – Cristal fotônicos otimizada pelo GA (a) para o modo TE12 e (b) para o Modo TM34. ........................................................................................................ 58 Figura 26 – Diagrama de banda de estrutura teste para o modo TE12. ..................... 59 Figura 27 – (a) Gráfico de evolução do ABC para o modo TE12 com 20 abelhas e 500 ciclos e (b) a célula unitária encontrada com PBG de 45,59%. .......................... 59 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107246 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107246 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107246 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107248 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107248 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107250 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107250 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107250 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107251 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107251 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107253 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107253 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107253 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107255 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107255 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107255 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107256 file:///C:/Users/FABIANO/Documents/backup_05-02-2020/mestrado/Disertação/Eviado%20Gilliard%20-%20Final/Dissertacao-Final.docx%23_Toc51107256 12 Figura 28 – Diagrama de banda de estrutura para o modo TE12. ............................. 60 Figura 29 – Diagrama de banda de estrutura para o modo TM34. ............................. 60 Figura 30 – Cristal fotônicas e célula unitária otimizada pelo ABC: (a) Modo TE12; (b) Modo TM34. ........................................................................................................ 60 Figura 31 – Gráfico de evolução do GA e ABC para o modo TE12. ........................... 61 Figura 32 – Gráfico de evolução do GA e ABC para o modo TM34. .......................... 62 Figura 33 – Diagrama de banda de estrutura para o modo TM12. ............................. 63 Figura 34 – Diagrama de banda de estrutura para o modo TE34. ............................. 63 Figura 35 – Cristais fotônicos e célula unitária otimizada pelo ABC: (a) Modo TM12; (b) Modo TE34..................................................................................................... 63 Figura 36 – Gráfico de evolução do ABC para o modo TM12. ................................... 64 Figura 37 – Gráfico de evolução do ABC para o modo TE34. .................................... 64 Figura 38 – Cristais fotônicos compondo as células unitárias encontrada na otimização com quantidade de abelhas totais igual a 20 para os seguintes modos: (a) TE12; (b) TE34; (c) TM12 e (d) TM34. .................................................. 65 Figura 39 – Gráficos de evolução para os quatros modos propostos com a quantidade de abelhas totais igual a 20. ............................................................ 66 Figura 40 – Gráfico comparativo entre a evolução do ABC e do AG para maximização da PBG do modo TE12. ................................................................. 66 Figura 41 – Gráfico comparativo entre a evolução do ABC e do AG para maximização da PBG do modo TM34. ................................................................ 67 13 LISTA DE TABELAS Tabela 1 – Parâmetros de controle do GA e ABC para otimização das funções bidimensionais. .................................................................................................. 47 Tabela 2 – Parâmetros de controle ABC para otimização do cristal fotônico. ........... 57 Tabela 3 – Parâmetros adotados para otimização com GA. ..................................... 58 Tabela 4 - Comparação de resultados entre GA e ABC ............................................ 62 Tabela 5 – Resultados obtidos na otimização com ABC com quantidade de abelhas totais igual 20 ..................................................................................................... 65 Tabela 6 – Resultados das PBGs obtidas por estruturas otimizadas pelo ABC. ....... 68 14 LISTA DE ALGORITMOS Algoritmo 1 – Pseudocódigo do GA clássico ............................................................. 24 Algoritmo 2 – Adaptação do ABC para o domínio discreto ....................................... 56 15 LISTA DE ABREVIATURAS E SIGLAS GA Genetic Algorithm ABC Artificial Bee Colony MBO Optimization with Marriage in Honey-Bees VBA Virtual Bee Algorithm BCO Bee Colony Optimization BSO Bees Swarm Optimization BA Bees Algorithm GABC Gbest-Guided Artificial Bee Colony PSO Particle Swarm Optimization MABC Modified Artificial Bee Colony CMOS Complementary Metal–Oxide–Semicondutor SOI Silicon-On-Insulator FEM Finite Element Method TE Transverse Electric TM Transverse Magnetic PBG Photonic Band Gap 16 LISTA DE SÍMBOLOS 𝒏 Índice de refração 𝜦 Período do arranjo cristalino 𝝀𝒐 Comprimento de onda no espaço livre �⃗⃗� Vetor campo elétrico �⃗⃗⃗� Vetor campo magnético �⃗⃗� Campo de deslocamento elétrico �⃗⃗� Vetor densidade de fluxo magnético 𝝆 Densidade de carga elétrica J Densidade de corrente �⃗� Vetor de posição no plano cartesiano 𝜺 Constante dielétrica relativa 𝜺𝟎 Constante dielétrica no espaço livre 𝝁 Permissibilidade magnética relativa 𝝁𝟎 Permissibilidade magnético no espaço livre 𝝎 Frequência angular 𝒄 Velocidade da luz no vácuo i Unidade imaginária �̂�, �̂�, �̂� Vetor unitário do plano cartesiano 17 SUMÁRIO 1 INTRODUÇÃO ......................................................................................... 20 2 META-HEURÍSTICAS BIOINSPIRADAS ................................................. 22 2.1 Algoritmo Genético (GA) ........................................................................ 22 2.2 Inteligência de Enxames ........................................................................ 25 2.3 Meta-Heurísticas Inspiradas nas Abelhas ............................................ 25 2.3.1 ARTIFICIAL BEE COLONY (ABC) ALGORITHM ..................................... 27 2.3.1.1 Adaptação do ABC para Otimizar de Problemas com Restrições ............ 30 2.3.1.2 Gbest-Guided ABC (GABC) ...................................................................... 31 2.3.1.3 Modified ABC (MABC) .............................................................................. 32 3 CRISTAIS FOTÔNICOS ........................................................................... 32 3.1 Equações de Maxwell ............................................................................. 34 3.2 Teorema de Bloch e Zona de Brillouin .................................................. 38 3.3 Banda Proibida Fotônica (Photonic Band Gaps - PBG) ...................... 41 3.4 Aplicação de Estruturas Cristalinas Fotônicas Bidimensionais ........ 42 3.4.1 GUIAS DE ONDAS ................................................................................... 43 3.4.2 FILTRO PASSA FAIXA ............................................................................. 45 3.4.3 DIVISOR DE POTÊNCIA .......................................................................... 45 4 METODOLOGIA ....................................................................................... 46 4.1 Comparação entre ABC e GA em Funções Multimodais Continuas. . 46 4.1.1 SIMULAÇÃO COMPARATIVA DA FUNÇÃO 1 ......................................... 48 4.1.2 ANÁLISE COMPARATIVA DA FUNÇÃO 2 ............................................... 50 4.2 Otimização de Estrutura Cristalina Fotônica Bidimensional .............. 54 4.3 Adaptação do ABC para Otimização do Cristal Fotônico ................... 54 5 RESULTADOS ......................................................................................... 57 5.1 Otimização com População Total de Abelhas Igual a 40 .................... 59 18 5.2 Otimização com População Total de Abelhas Igual a 20 .................... 64 6 CONSIDERAÇÕES FINAIS ..................................................................... 68 19 Publicações associadas à dissertação MALHEIROS-SILVEIRA, G. N.; DELALIBERA, F. G. Inverse design of photonic structures using an artificial bee colony algorithm. Applied Optics, v. 59, n. 13, p. 4171-4175, 2020. MALHEIROS-SILVEIRA, G. N.; DELALIBERA, F. G.; RODRIGUEZ-ESQUERRE, V. F. Photonic bandgap optimized by artificial bee colony algorithm. In: Nanoengineer- ing: Fabrication, Properties, Optics, Thin Films, and Devices XVII. International Soci- ety for Optics and Photonics, 2020. p. 1146709. 20 1 INTRODUÇÃO Problemas de otimização estão presentes em diversas áreas do cotidiano tais como: economia, engenharia, processos industriais, entre outros. Um bom exemplo pode ser expresso quando em um processo industrial, no qual, mediante a variação de determinados fatores ou variáveis, busca-se a otimização de um determinado produto pela redução do seu custo/perdas ou o aumento da sua produção. As ferramentas de otimização são empregadas a fim de otimizar tais problemas de modo computacional. Estas ferramentas utilizam uma função objetivo, podendo ser contínua ou discreta, para quantificar as possíveis soluções do problema. Podemos dividir essas ferramentas em dois grupos de programação: linear e não-linear (LUENBERGER, 1984). A programação linear é utilizada na otimização de problemas cuja a função objetivo é composta por equações ou inequações lineares, de 1º grau. Um exemplo deste método é o Simplex, apresentando por Dantzig, em seu trabalho publicado em 1953 (DANTZIG, 1953). Com o intuito de proporcionar uma melhor compreensão, a programação não- linear pode ser subdividida em três métodos, sendo eles: determinísticos, enumerativos e estocásticos. Estes métodos são aplicados a problemas mais complexos, geralmente representados por equações diferenciais (LUENBERGER, 1984). Os métodos determinísticos são baseados nas informações decorrentes do gradiente da função objetivo, buscando os pontos onde se anula ou muda de sinal. A utilização deste método é limitada pela complexidade da função objetivo, podendo ser aplicável somente às funções contínuas e deriváveis em todo o seu domínio. Os métodos enumerativos têm a capacidade de encontrar a melhor solução para um problema, porém sua dinâmica de busca consiste em uma varredura completa de todas as suas possíveis soluções, tornando-se extremamente exaustivo. A implementação deste método é limitada pela quantidade de possíveis soluções e recursos computacionais. Um exemplo do método enumerativo foi proposto para otimização de funções discretas por (MARKOWITZ, 1957). Os métodos estocásticos ou heurísticos são algoritmos probabilístico que fazem a busca pelo ótimo de modo aleatório, ordenado pelos seus parâmetros de 21 controle. Estes métodos não garantem um ótimo absoluto, porém buscam uma aproximação eficiente com recursos computacionais e tempo reduzidos. Métodos heurísticos mais sofisticados vêm sendo cada vez mais implementados, sob a nomenclatura de meta-heurísticas, nas quais uma mesma heurística pode ser aplicada de modo genérico a diferentes problemas de otimização (LUKE, 2013). A maioria dos problemas a serem otimizados, no meio cientifico, são reconhecidamente complexos. Desse modo, surge então a necessidade de utilizar as meta-heurísticas. Tratando-se de eletromagnetismo, encontram-se diversos problemas com elevados níveis de complexibilidade, dentre os quais encontra-se a otimização de estruturas cristalinas fotônicas, chamadas de cristais fotônicos, que é o tema sobre o qual desenvolve-se este trabalho. Cristais fotônicos são arranjos periódicos de diferentes materiais dielétricos que exibem algumas propriedades referentes a propagação de ondas eletromagnéticas. Dentre as suas propriedades está a banda fotônica proibida (PBG), na qual o cristal, à uma determinada faixa de frequência, atua como um espelho “perfeito”, bloqueando a propagação de luz. A PBG pode ser maximizada alterando o arranjo, os materiais e/ou a geometria das células unitárias de um cristal fotônico. Inúmeros algoritmos foram aplicados para otimizar os cristais fotônicos, sendo alguns deles: • Algoritmo Evolutivo (EA) (PREBLE; LIPSON; LIPSON, 2005); • Algoritmos Genético (GA) (HAKANSSON; SÁNCHEZ-DEHESA; SANCHIS, 2005) e (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE; HERNANDEZ-FIGUEROA, 2011); • Otimização por Enxame de Partículas (PSO) (DJAVID; MIRTAHERI; ABRISHAMIAN, 2009); • Sistemas Imunológicos Artificiais (AIS) (FERREIRA; MALHEIROS- SILVEIRA; HERNANDEZ-FIGUEROA, 2015) e (SISNANDO et al., 2015); • Evolução Diferencial (ED) (BOR; TURDUEV; KURT, 2016); • Otimizador Gray Wolf (GWO) (MIRJALILI et al., 2017). Este trabalho propõe a otimização de estruturas cristalinas fotônicas bidimensionais utilizando uma nova meta-heurística, apresentada há menos de duas 22 décadas, por Dervis Karaboga, denominada: Algoritmo Colônia de Abelhas Artificiais (Algorithm Bee Colony – ABC) (KARABOGA, 2005). O objetivo é maximizar a PBG de estruturas compostas de Silício (Si) e ar, em um arranjo triangular, alterando somente a sua geometria. Para a validação dos resultados serão utilizados valores de PBG de estruturas otimizadas por outra meta-heurística bastante conceituada, o GA (MALHEIROS- SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007). 2 META-HEURÍSTICAS BIOINSPIRADAS Meta-heurísticas bioinspiradas são algoritmos inspirados em fenômenos da natureza. Dentre a diversidade de meta-heurísticas bioinspiradas existentes na literatura, este capítulo faz destaque somente a duas: Algoritmo Genético (GA) e Algoritmo Colônia de Abelhas Artificiais (ABC), dos quais será introduzido uma breve revisão a seguir. 2.1 Algoritmo Genético (GA) O GA é um método de busca estocástica inspirado na teoria da origem das espécies através da seleção natural (Darwin, 1859). Foi proposto por John Holland e seus colegas da Universidade de Michigan, nos anos de 1960, e publicando em seu livro Adaptation in Natural and Artificial Systems (HOLLAND, 1975). Outros algoritmos evolutivos baseados na teoria da evolução Darwiniana foram apresentados nesse período, como de Fogel, Owens e Walsh (FOGEL, 1966). Entretanto, o GA tem se destacando cada vez mais no meio cientifico, tendo grande relevância até os dias de hoje. Classificado como uma meta-heurística populacional, o GA é representado por um número finito de cromossomos (indivíduos), sendo cada indivíduo uma possível solução do problema. A interação do algoritmo ocorre pelo cruzamento (crossover) dos cromossomos, gerando uma nova geração de indivíduos, fazendo uma analogia ao cruzamento genético proposto nas leis de Mendel. Cada cromossomo é representado por uma quantidade de genes pré- estabelecida, representados computacionalmente pelas variáveis do problema (MITCHELL, 1998). Inicialmente, a população é gerada de modo aleatório, para só então avaliar a 23 aptidão de cada indivíduo através de uma função fitness ou função custo. Através dos valores apresentados, cada indivíduo apresenta um percentual de aptidão, podendo ser calculado pela equação (1). Quanto maior o percentual, maior a probabilidade do indivíduo ser escolhido para gerar novos indivíduos, ocorrendo então a seleção natural apresentada por Darwin. 𝑃𝑖 = 𝐹𝑖 ∑ 𝐹𝑛 𝑁 𝑛=1 . 100% (1) Onde Pi é a probabilidade de seleção e Fi é o valor da função objetivo do indivíduo de índice i. N é o número total de indivíduos em cada geração. Um exemplo do cruzamento entre dois indivíduos é ilustrado na Figura 1, onde dois cromossomos denominados “pais” são divididos em duas ou mais partes entre o(s) ponto(s) de crossover e a combinação das partes de cada cromossomo gera novos cromossomos, denominados “filhos”. Além do cruzamento, um ou mais genes do cromossomo filho pode sofrer alterações pelo operador mutação, alterando os seus valores de modo aleatório, simulando a modificação genética do indivíduo. Um mesmo problema pode ser otimizado pelo GA, obtendo diferentes resultados, pois além de ser um algoritmo probabilístico, ele depende de alguns parâmetros de controle ajustáveis que interferem diretamente no resultado, a depender do problema, como: • Tamanho da população: Quantidade de indivíduos que compõe a população a cada geração, podendo ser estático ou dinâmico, ao decorrer de cada Figura 1 – Exemplo de operação crossover entre dois indivíduos, com apenas um ponto de crossover e um bit de mutação. 24 geração. Um número baixo de indivíduos ocasiona uma limitação na exploração do espaço de busca, enquanto que uma quantidade muito grande exige um tempo computacional elevado. • Taxa de crossover (cruzamento): Define a porcentagem entre os cromossomos que serão sujeitos ao cruzamento e os que serão preservados para próxima geração, podendo variar entre 0% a 100%. Sendo que, quando a probabilidade é de 100%, todos os cromossomos são sujeitos ao cruzamento, e quando é 0%, toda a população é replicada para próxima geração podendo ser alteradas somente pela mutação. Normalmente é estabelecido valores entre 85% a 95%, permitindo que entre 15% a 5% dos indivíduos mais aptos façam parte da geração seguinte. • Taxa de mutação: É a porcentagem de genes de um cromossomo sujeita a mutação, sendo os genes escolhidos de modo aleatório. Esta taxa pode variar de 0% (quando nenhum gene é sujeito a mutação) à 100% (quando todos os genes são alterados, ocorrendo a inversão de um cromossomo). Normalmente, são estabelecidos valores baixos entre 0,5% a 1%. • Método de Seleção: A seleção dos cromossomos “pais” pode ser obtida através de diversos métodos, como, por exemplo: Roleta, Elitismo, Seleção de Boltzmann, Seleção de Torneiro, entre outros. A análise detalhada de alguns dos métodos de seleção pode ser encontrada em (MITCHELL pág. 124, 1998). Este trabalho abordara apenas o método de roleta, no qual a seleção dos indivíduos “pais” ocorre através de um sorteio aleatório considerando cada indivíduo da geração atual, tendo maior probabilidade de serem sorteados os indivíduos com maior aptidão. A porcentagem de cada indivíduo é calculada por meio à equação (1). Tendo em vista o exposto até então, o pseudocódigo do GA clássico pode ser resumido como: Algoritmo 1 – Pseudocódigo do GA clássico. 1 - Geração de População inicial 2 - Avaliação 3 - 4 - Enquanto Condição de parada Fazer: 5 - Seleção 6 - Crossover 7 - Mutação 8 - Avaliação 9 - Fim 25 2.2 Inteligência de Enxames Algoritmos de Inteligência de Enxames são inspirados no comportamento do trabalho coletivo entre grupos de animais ou insetos, como: pássaros, peixes, abelhas, formigas, entre outros. Tais algoritmos podem ser fundamentados por dois conceitos principais: Auto-organização e divisão de trabalho (BONABEAU, 1999). • A auto-organização corresponde ao conjunto de mecanismos dinâmicos que institui, de modo hierárquico, todo um sistema e é formado por regras básicas de relacionamento entre cada parte. São apontadas por Bonabeau (1999) quatro propriedades básicas para troca de informação: análise quantitativa sobre as fontes de alimentos encontradas (feedback positivo); informações de saturação de alimento em determinadas fontes (feedback negativo); busca aleatória por novas fontes de alimentos (flutuação) e quantidade mínima de indivíduos para interações aleatórias (múltiplas interações). • Divisão de trabalho é o conceito que estabelece que cada indivíduo é direcionado para a realização de tão somente uma tarefa, tornando-se um indivíduo especializado, de modo a permitir que todas as tarefas sejam realizadas simultanea- mente (paralelismo), tornando o sistema mais eficiente do que se um mesmo indivíduo efetuar todas as tarefas sequencialmente. O conceito de paralelismo reduz o tempo de produção e energia dispendida pelos indivíduos, uma vez que não exige que um mesmo indivíduo troque de função, permitindo que cada um se responsabiliza so- mente por uma tarefa (BONABEAU, 1999). 2.3 Meta-Heurísticas Inspiradas nas Abelhas Estudos a respeito das características naturais de abelhas produtoras de mel inspirou o desenvolvimento de alguns algoritmos de otimização. O comportamento de tomada de decisão coletiva na busca por melhores fontes de alimentos (SEELEY, 1991) e métodos de acasalamento das abelhas rainhas (ADAMS et al., 1977) foram a base para estes algoritmos. Inicialmente, algoritmos inspirados em abelhas foram baseados no comporta- mento de procriação das abelhas rainhas para gerarem novas abelhas (linhagem). Esses algoritmos utilizam métodos probabilísticos para simular a procriação das abe- lhas, tais como: 26 • Otimização com Casamento em Abelhas (MBO – Optimization with Marriage in Honey-Bees) (ABBASS, 2001); • Otimização de Acasalamento das Abelhas (HBMO – Honey-Bees Mating Opti- mization) (HADDAD, 2006). Posteriormente, foram apresentados algoritmos baseados no comportamento de busca por melhores fontes de alimentos (forrageamento), dentre os quais: • Colmeia (BeeHive) (WEDDE et al., 2004); • Algoritmo Colônia de Abelhas (ABC - Artificial Bee Colony) (KARABOGA, 2005); • Algoritmo de Abelhas Virtuais (VBA - Virtual Bee Algorithm) (YANG, 2005); • Otimização de Colônias de Abelhas (BCO - Bee Colony Optimization) (TEODOROVIC, 2005); • Otimização de Enxame de Abelhas (BSO - Bees Swarm Optimization) (DRIAS et al., 2005); • Algoritmo de Abelhas (BA - Bees Algorithm) (PHAM et al., 2006). Uma comparação relacionada a quantidade de publicações referente aos algo- ritmos HBMO, BCO, BA e ABC (realizada no período entre 2005 a 2011) apontou que 54% dos pesquisadores demonstraram preferência pela abordagem utilizando o ABC em comparação as demais, conforme ilustra a Figura 2. Isso indica que o ABC vem se destacando cada vez mais entre as meta-heurísticas inspiradas no comportamento das abelhas (KARABOGA, 2014). Fonte: Adaptado de KARABOGA (2014). Dentre as publicações citadas na Figura 2, o ABC foi a opção para resolver problemas em diversas áreas. Em lista apresentada por Karaboga (2014), foram doze as áreas nas quais mais se apurou publicações relacionadas a utilização do ABC, 54% 15% 18% 13% ABC HBMO BA BCO Figura 2 – Porcentagem de publicação relacionada a cada Algoritmo HBMO, BA, BCO e ABC. 27 durante os cinco primeiros anos logo após a sua publicação, conforme a distribuição gráfica demonstrada na Figura 3. Fonte: (KARABOGA, 2014). 2.3.1 ARTIFICIAL BEE COLONY (ABC) ALGORITHM O ABC foi proposto como um novo algoritmo de inteligência de enxames para otimização de funções multimodais contínuas, buscando o mínimo global de três funções previamente conhecidas. Karaboga apresentou o ABC como um algoritmo mais simples e flexível em comparação aos demais algoritmos de enxames (KARABOGA, 2005). No ABC, o conceito de divisão de trabalho ocorre classificando as abelhas em três diferentes grupos, sendo eles: Abelhas empregadas (employed), abelhas espectadoras (onlooker) e as abelhas exploradoras (scouts). As abelhas empregadas são enviadas às fontes de alimentos e fornecem os feedbacks positivos e negativos à colmeia. As espectadoras fazem a busca por novas fontes ao redor das abelhas empregadas, mediante as informações de seus feedbacks, e as exploradoras são enviadas a realizar uma busca por novas fontes de alimento, de modo aleatório, para aumentar a diversidade de busca. A posição de cada abelha é relacionada a uma fonte de alimento, aonde a quantidade de néctar em cada fonte é avaliada quantitativamente mediante a uma função objetivo (fitness). Cada possível solução do problema representa uma fonte de alimento. 7% 22% 3% 13% 15% 4% 5% 6% 10% 10% 3% 2% Redes neurais Engenharia Industrial Engenharia Mecânica Engenharia Elétrica Engenharia Eletrônica Engenharia de Controle Engenharia Civil Engenharia de Software Processamento de Imagens Mineração de dados Redes de sensores Estrutura proteica Figura 3 – Publicações em diferentes áreas nos cinco primeiros anos do algoritmo ABC. 28 Um dos motivos do ABC ser apresentado como simples pode se dar pelo fácil manuseio dos parâmetros de controle. No ABC clássico são apresentados quatro parâmetros: Quantidade de abelhas empregadas, quantidade de abelhas espectadoras, condição de parada e o parâmetro limite. O parâmetro limite define a quantidade máxima de ciclos na qual uma abelha empregada pode ficar em uma mesma fonte de alimentos, sem que seu alimento se esgote. Atingido esse limite, a abelha empregada estará sujeita a se tornar uma abelha exploradora, até que encontre outra fonte de alimento (KARABOGA, 2005). Para facilitar a compreensão do ABC, uma síntese de seu comportamento de busca pode ser representada em três etapas, conforme ilustra o fluxograma na Figura 4. • Etapa 1 – Cada abelha empregada é enviada, aleatoriamente, para uma fonte de alimento e calcula a quantidade de néctar através da função objetivo. A informação de cada fonte de alimento é compartilhada com a colmeia. Esta etapa ocorre somente uma vez, para inicialização do algoritmo. • Etapa 2 – Cada abelha espectadora escolhe uma abelha empregada, semelhante ao método de roleta, para realizar a busca por uma nova fonte de alimentos ao seu derredor, por meio da chamada “dança das abelhas”. Se a fonte de alimento encontrada obtiver a quantidade néctar maior do que a da abelha empregada escolhida, então a abelha empregada abandona a sua fonte de alimento e se move para a nova fonte de alimento encontrada pela abelha espectadora. • Etapa 3 – A cada ciclo, uma abelha empregada cuja a fonte de alimento se esgotou, torna-se uma abelha exploradora, sendo enviada a uma nova fonte de alimentos, aleatoriamente, gerando uma diversidade entre o campo de busca. Ao final de cada ciclo, a posição da abelha com a melhor fonte de alimento encontrada é salva em uma variável, usualmente chamada de Gbest, e uma abelha empregada com fontes de alimento esgotada é separada para se tornar uma exploradora no próximo ciclo. Repete-se as etapas 2 e 3 a cada ciclo, sendo a quantidade de ciclo dada pela condição de parada, que são usualmente determinadas pela quantidade de ciclos, tempo ou valor a ser encontrado. 29 Figura 4 – Fluxograma do comportamento do ABC separado em etapas. As abelhas reais fazem a busca ao redor da fonte de alimento através do movimento chamado “dança das abelhas”. Para representar este movimento, as abelhas virtuais realizam a busca por meio da equação abaixo: 𝑦𝑖𝑗 = 𝑥𝑖𝑗 + Ө(𝑥𝑖𝑗 − 𝑥𝑘𝑗) (2) 30 Onde yij e xij representam as abelhas espectadoras e empregadas, respectivamente; i indica o índice de cada abelha e j cada variável do problema (dimensão); xkj é a abelha empregada de índice k e dimensão j, onde k é escolhido aleatoriamente e diferente de i. θ é apresentado por Karaboga como um parâmetro de aceleração, normalmente um número aleatório entre [-1,1]. A probabilidade de envio de abelhas espectadora a cada fonte de alimentos pode ser representada, probabilisticamente, através da equação (1), enviando um número maior de abelha às fontes de alimentos com maior quantidade de néctar. Onde Fi é o valor da função objetivo da abelha de índice i e N é o número total de abelhas empregadas. Submetido a uma análise comparativa com outras meta-heurísticas, o ABC mostrou ser eficiente para resoluções de problemas de engenharia, tanto para funções multimodais com múltiplas variáveis (KARABOGA, 2007), como para funções multidimensionais (KARABOGA, 2008). O ABC tem recebido atenção e interesse crescentes por parte de alguns pesquisadores, que têm proposto algumas variações na versão clássica apresentada por Karaboga, com objetivo de melhorar a sua eficiência e ampliar sua aplicação. Algumas dessas proposições serão listadas a seguir. 2.3.1.1 Adaptação do ABC para Otimizar de Problemas com Restrições O ABC clássico apresentou limitações por ser aplicar somente a otimização de funções contínuas. Karaboga e Basturk, com o objetivo de ampliar sua aplicação à problemas com restrições no domínio (KARABOGA; BASTURK, 2007), propuseram algumas adaptações utilizando o método de seleção apresentado por Deb (2000) e uma modificação na geração de novas soluções do problema. Basicamente, o método de Deb consiste em um operador de seleção por tor- neio, gerando duas soluções, simultaneamente, e comparando-as, obedecendo a três critérios: • Soluções viáveis tem preferência se comparadas a soluções inviáveis; • Entre duas soluções viáveis, qual tiver o maior valor da função objetivo tem preferência; 31 • Entre duas soluções inviáveis, tem preferência a com menor violação de restrição. A modificação na busca por novas fontes de alimento se deu mediante a uma alteração da equação (2), propondo uma condição de comparação entre fatores para considerar também alguns valores da própria abelha empregada, conforme mostra a equação (3). 𝑦𝑖𝑗 = { 𝑥𝑖𝑗 + Ө(𝑥𝑖𝑗 − 𝑥𝑘𝑗), 𝑆𝑒 𝑅𝑗 < 𝑀𝑅 𝑥𝑖𝑗 , 𝑆𝑒𝑛ã𝑜 (3) Onde Rj é um número aleatório entre [0,1] para cada dimensão j e MR é cha- mado de taxa de modificação, um parâmetro de controle que determina se xij será modificado ou não, representado por valores entre [0,1]. Com isso, é possível modifi- car, de modo aleatório, somente as variáveis em algumas dimensões da abelha em- pregada, gerando maiores configurações possíveis para a nova solução. 2.3.1.2 Gbest-Guided ABC (GABC) Zhu e Kwong relataram que o ABC é muito eficiente para buscas aleatórias, porém apresenta limitações na exploração intensificada em uma busca local (ZHU; KWONG, 2010). Propuseram então uma alteração na equação (2) do ABC clássico, baseada no algoritmo de Otimização por Exame de Partícula (PSO – Particle Swarm Optimization) (EBERHART; KENNEDY, 1995), acrescentando os valores referente a melhor solução encontrada dada pela variável Gbest, conforme mostra a equação (4). 𝑦𝑖𝑗 = 𝑥𝑖𝑗 + Ө(𝑥𝑖𝑗 − 𝑥𝑘𝑗) + 𝛹𝑖𝑗(𝑔𝑗 − 𝑥𝑖𝑗) (4) Onde gj é o valor da melhor solução na dimensão j e Ψij um número aleatório entre [0, C], sendo C uma constante não negativa. O experimento de comparação mostrou melhores resultados do GABC em comparação ao ABC para seis funções multimodais e multidimensionais contínuas (ZHU; KWONG, 2010). 32 2.3.1.3 Modified ABC (MABC) Uma nova proposta pra aumentar a performance do ABC, apresentou modifi- cações nos métodos de: busca por fontes vizinhas (2), seleção das abelhas especta- doras (1) e busca aleatória das abelhas exploradoras (HADIDI et al., 2010). O novo método de busca por fontes vizinhas propõe uma escolha aleatória do uso da equação (2), na qual uma constante dada por um parâmetro Pa com valores entre [0,1], é comparada a um número aleatório entre [0,1] para cada parâmetro j de uma abelha. Se o número aleatório for menor que Pa, a nova solução na dimensão j é dada pela equação (2), senão atribui o valor do parâmetro j da abelha empregada. O método de seleção das abelhas espectadoras propõem a seleção por torneio. Um número M de fontes de alimentos são selecionadas aleatoriamente das abelhas empregadas, então a fonte escolhida pela abelha espectadora será a com maior quan- tidade de néctar. Desta maneira, não se torna mais necessário utilizar a equação (1). Também relataram que a busca aleatória das abelhas exploradoras se torna eficiente somente na inicialização da busca. Propuseram, portanto, um método deno- minado mais eficiente para busca aleatória, considerando o desvio padrão de cada abelha conforme apresenta na equação (5). 𝑦𝑖𝑗 = 𝑥𝑖𝑗 + 𝑁(0, 𝜎𝑗) (5) Onde, yij é a nova fonte encontrada pela abelha espectadora, σj é o desvio pa- drão para cada uma das variáveis da dimensão j e N(0, σj) é um número aleatório, normalmente distribuído com média 0 e desvio padrão de σj. Gao e Liu utilizaram vinte e seis funções benchmark escaláveis com até 100 dimensões (j = 100), e outras duas com até 300 dimensões (j = 300) para validação do algoritmo, comparando-o com o ABC clássico. O MABC apresentou melhores re- sultados na maioria das simulações comparadas ao ABC (GAO; LIU, 2012). 3 CRISTAIS FOTÔNICOS Cristais fotônicos (Photonic Crystals) são chamados de cristais por serem estruturas periódicas cristalinas e fotônicos por afetarem os movimentos dos fótons 33 (luz). Esta estrutura é composta de materiais com índices de refração distintos, em arranjos quadrados ou triangulares (JOANNOPOULOS et al., 1995). Os cristais fotônicos são classificados como unidimensional (1-D), bidimensio- nal (2-D) e tridimensional (3-D). Estruturas com periodicidade em apenas uma direção são denominadas unidimensionais, já aquelas com periodicidade em um plano são denominadas bidimensionais. Quando a periodicidade se estende ao longo dos três eixos cartesianos (x, y, z), tem-se os chamados cristais fotônicos tridimensionais, con- forme ilustrado na Figura 5. Fonte: (JOANNOPOULOS et al., 1995. p. 4). Um exemplo de estrutura cristalina fotônica unidimensional é apresentado na Figura 6 (a), formada por um arranjo periódico de dois materiais dielétricos com periodicidade Λ relacionado a uma constante de rede a na direção �̂�, com altura H e índices de refração n1 e n2, suspenso numa estrutura com índice de refração n3 onde normalmente: n1 > n2, n3. A relação entre o período Λ e o comprimento de onda no espaço livre (λ0), permite que o cristal se comporte em três principais regimes: • Sub-comprimento de ondas (sub-wavelength) (λ0 ≫ Λ) - a estrutura se comporta como um guia de onda, permitindo a propagação da luz (HALIR et al., 2015). • Reflexão (λ0 / 2nef = Λ) - a estrutura se comporta como um espelho “perfeito” (Bragg mirror), onde a luz é refletida de volta ao emissor; • Difração (λ0 ≪ Λ) - a luz é irradiada em diferentes direções ao redor da estrutura cristalina, ocasionando a dispersão no meio; A propagação do campo elétrico nos três regimes apresentados acima é ilustrada na Figura 6 (b). Figura 5 – Exemplo de cristais fotônicos unidimensional, bidimensional e tridimensional. 34 Figura 6 – (a) Estrutura cristalina com periodicidade discreta Λ na direção �̂� e (b) propagação do campo �⃗� (x,z) entre a estruturas para os três principais regimes: sub-comprimento de ondas, reflexão e difração. Fonte: (HALIR et al., 2015). Mediante estes regimes, os cristais fotônicos vêm sendo estudados para uma infinidade de aplicações em dispositivos ópticos integrados, como em sensores (ZHANG; ZGAO; Ri-Ging, 2015), laser, guias de ondas, acopladores, entre outros (DUTTA et al., 2016). Com as crescentes evoluções nas telecomunicações, por meio da transmissão de dados por sinais luminosos com comprimentos de onda entre 1,3μm a 1,6μm, os cristais fotônicos têm sido projetados para atuarem na região de infravermelho, tornando necessária a fabricação de nano cristais fotônicos com período Λ menor ou igual a 0,5μm. Além disso, estudos de fabricação utilizando a mesma tecnologia utilizada para fabricação de microeletrônica, complementary metal–oxide– semicondutor (CMOS) em wafers de silício (silicon-on-insulator - SOI), para nano cristais fotônicos têm apresentado resultados satisfatórios para cristais fotônicos bidimensionais (BOGAERTS et al., 2005). Neste capítulo, será apresentada uma análise básica baseada nas equações de Maxwell do comportamento dos campos elétricos e magnéticos e uma análise de estrutura de banda proibida em cristais fotônicos unidimensionais. 3.1 Equações de Maxwell A propagação de luz, bem como todos fenômenos eletromagnéticos, são regidos pelas quatro equações macroscópicas de Maxwell, sendo elas: (a) (b) 35 ∇ ∙ �⃗� = 0 ∇ × �⃗� + 𝜕�⃗� 𝜕𝑡 = 0 (6) ∇ ∙ �⃗⃗� = 𝜌 ∇ × �⃗⃗� − 𝜕�⃗⃗� 𝜕𝑡 = 𝐽 Onde �⃗� e �⃗⃗� são os vetores campo elétrico e campo magnético, respectivamente, �⃗⃗� é o campo de deslocamento elétrico, �⃗� é o vetor densidade de fluxo magnético, ρ é a densidade de carga e J a densidade de corrente, sendo todos dependente do tempo e espaço. Restringido a propagação de luz em um meio dielétrico homogêneo e sem perdas, pertencentes a um plano cartesiano com vetor de posição 𝑟 , nos permite fazer algumas aproximações para simplificar a análise: • O meio é ausente de fontes de luz internas, podendo assumir que a densidade de carga e densidade de corrente são nulas (ρ = J = 0); • Os campos são “fracos”, podendo então considerar somente os seus regimes lineares; • O material é macroscópico e isotrópico, portanto, �⃗� e �⃗⃗� podem ser relacionadas pela constante dielétrica no espaço livre ɛ0 e uma função dielétrica escalar ɛ(𝑟 ,t), também chamada de permissibilidade relativa. Logo: �⃗⃗� = ɛ0.ɛ(𝑟 ,t).�⃗� ; • A dependência explicita da constante dielétrica em relação ao tempo é nula ɛ(𝑟 ,t) = ɛ(𝑟 ). • A constante de permissibilidade magnética relativa μ é igual a unidade, portanto �⃗� = μ0.μ(𝑟 ).�⃗⃗� (𝑟 ) = μ0.�⃗⃗� (𝑟 ). onde, μ0 = 4π.10-7 Henry/m. É pertinente dizer que, mediante essas aproximações, são desconsiderados alguns fenômenos físicos dos materiais. No entanto, muitas propriedades úteis e elementares dos materiais surgem do caso elementar de materiais lineares e sem perdas. Deste modo, a teoria de cristais fotônicos fica mais simples, além de proporcionar uma base às teorias mais sofisticadas (JOANNOPOULOS et al., 1995). Deste modo, as quatro equações de Maxwell podem ser reescritas como: ∇ ∙ �⃗⃗� (𝑟 , 𝑡) = 0 ∇ × �⃗� (𝑟 , 𝑡) + 𝜇0 𝜕�⃗⃗� (𝑟 , 𝑡) 𝜕𝑡 = 0 (7) ∇ ∙ [𝜀(𝑟 )�⃗� (𝑟 , 𝑡)] = 0 ∇ × �⃗⃗� (𝑟 , 𝑡) − 𝜀0𝜀(𝑟 ) 𝜕�⃗� (𝑟 , 𝑡) 𝜕𝑡 = 0 36 O campo elétrico e magnético são funções complicadas tanto no domínio do tempo como do espaço. Sabendo que as equações de Maxwell são lineares, torna-se possível separar a dependência temporal da espacial, por meio de modos harmônicos da seguinte forma: �⃗⃗� (𝑟 , 𝑡) = �⃗⃗� (𝑟 )𝑒−𝑖𝜔𝑡 (8) �⃗� (𝑟 , 𝑡) = �⃗� (𝑟 )𝑒−𝑖𝜔𝑡 (9) Onde os campos variam de modo senoidal em relação ao tempo com frequência ω e o campo físico é obtido extraindo a parte real da equação. Substituindo �⃗⃗� (𝑟 ,t) e �⃗� (𝑟 ,t) nas quatro equações de Maxwell, obtemos: ∇ ∙ �⃗⃗� (𝑟 ) = 0 (10) ∇ ∙ [𝜀(𝑟 )�⃗� (𝑟 )] = 0 (11) ∇ × �⃗� (𝑟 , 𝑡) + 𝜇0 𝜕�⃗⃗� (𝑟 )𝑒−𝑖𝜔𝑡 𝜕𝑡 = 0 (12) ∇ × �⃗⃗� (𝑟 , 𝑡) − 𝜀0𝜀(𝑟 ) 𝜕�⃗� (𝑟 )𝑒−𝑖𝜔𝑡 𝜕𝑡 = 0 (13) As equações (10) e (11) são condições de transversalidade, nas quais uma onda plana dada por �⃗⃗� (𝑟 ) = a. 𝑒𝑖�⃗� .𝑟 para um dado vetor de onda �⃗� , será fisicamente aceitável somente se a.�⃗� = 0. Derivando no tempo o campo magnético �⃗⃗� (𝑟 )𝑒−𝑖𝜔𝑡 e o campo elétrico �⃗� (𝑟 )𝑒−𝑖𝜔𝑡 das equações (12) e (13), obtemos: ∇ × �⃗� (𝑟 ) − 𝑖𝜔𝜇0�⃗⃗� (𝑟 ) = 0 (14) ∇ × �⃗⃗� (𝑟 ) − 𝑖𝜔𝜀0𝜀(𝑟 )�⃗� (𝑟 ) = 0 (15) 37 Dividindo a equação (15) pela permissibilidade dielétrica ɛ(𝑟 ) e tomando o rotacional, torna-se possível utilizar a equação (14) para eliminar o parâmetro campo elétrico �⃗� (𝑟 ), obtendo uma equação com somente o campo magnético �⃗⃗� (𝑟 ) apresentada na equação (16). As constantes ɛ0 e 𝜇0 podem ser relacionadas com a velocidade da luz no vácuo c = 1 √𝜀0𝜇0 . ∇ × ( 1 𝜀(𝑟 ) ∇ × �⃗⃗� (𝑟 )) = ( 𝜔 𝑐 ) 2 �⃗⃗� (𝑟 ) (16) A equação (16) é chamada de Equação Mestre. Sendo sujeita às condições de transversalidade, apresenta todas informações sobre o campo magnético �⃗⃗� (𝑟 ). Podendo, para uma dada estrutura com uma função dielétrica escalar ɛ(𝑟 ), ser resolvida no domínio da frequência, encontrando as componentes harmônicas do campo magnético, também chamadas de modos ou estado. Para obter informações sobre o campo elétrico, basta substituir o valor de �⃗⃗� (𝑟 ) na equação (15), obtendo a seguinte equação: �⃗� (𝑟 ) = 𝑖 𝜔𝜀0𝜀(𝑟 ) ∇ × �⃗⃗� (𝑟 ) (17) Também é possível obter o campo elétrico por meio de uma única equação, porém torna-se mais complexa sua resolução pelo fato de os componentes não tangenciais a uma interface dielétrica serem descontínuos, motivo pelo qual a equação (16) é usualmente utilizada. A Equação Mestre pode ser tratada como um autovetor, onde o resultado do autovetor �⃗⃗� (𝑟 ) multiplicado por um operador Θ̂ é dado pela própria função multiplicada por uma constante (ω.c)2. Esta constante é chamada de autovalor. Θ̂�⃗⃗� (𝑟 ) = ( 𝜔 𝑐 ) 2 �⃗⃗� (𝑟 ) (18) Θ̂ = ∇ × ( 1 𝜀(𝑟 ) ∇ ×) (19) 38 Onde o operador Θ̂ é linear, simétrico e também real e positivo, pois todas constantes dielétricas ɛ(𝑟 ) são reais e positivas nos comprimentos de onda utilizados na óptica. Portanto, o operador contém todas as informações sobre a distribuição espacial da constante dielétrica, podendo ser associado ao operador Hermitiano, onde seu valor é igual ao seu transposto conjugado (Θ̂ = Θ̂ ∗ ). Outra importante informação do eletromagnetismo em meios dielétrico é que não possui uma escala de comprimento fundamental, tornando a Equação Mestre re- escalável. Portanto, dada uma estrutura com constante dielétrica ɛ(𝑟 ) a uma determinada frequência ω, obtém-se os mesmos modos harmônicos �⃗⃗� para outra estrutura com constante dielétrica ɛ(𝑟 )/s2 se as frequências forem re-escaladas para (ω.s). Isto permite que uma estrutura microscópica possa ser analisada por meio de outra, com escala maior. 3.2 Teorema de Bloch e Zona de Brillouin Uma descrição conveniente do comportamento da luz em estruturas periódicas pode ser expressa mediante uma analogia com a física do estado sólido. Felix Bloch descreveu que o movimento dos elétrons em um sólido pode ser expresso como um produto de uma onda plana de um elétron livre com uma função periódica, chamando então de teorema de Bloch (KITTEL, 2004). Em cristais fotônicos, o sistema translacional ao longo da periocidade dos cristais apresenta uma relação simétrica da constante dielétrica no sentido da periodicidade, permitindo expressar pela seguinte equação: 𝑇�̂�𝜀(𝑟 ) = 𝜀(𝑟 ) = 𝜀(𝑟 + 𝑑 ) (20) Podendo ser relacionado por um conjunto de vetores múltiplos inteiro da distância periódica: 𝜀(𝑟 ) = 𝜀(𝑟 + �⃗� ) (21) Onde os vetores �⃗� são chamados de rede de Bravais ou rede no espaço real. Por exemplo, o conjunto de vetores da rede de Bravais de uma estrutura 39 unidimensional ilustrada na Figura 6 (a) é dado por: �⃗� = γ.Λ. �̂�. Logo, para uma estrutura tridimensional, o conjunto de vetores é dado por: �⃗� = α.Λ. �̂� + β.Λ. �̂� + γ.Λ. �̂�, onde α, β e γ são números inteiros. γ Além da rede no espaço real, uma estrutura periódica é formada por um conjunto de vetores de ondas planas representado geralmente pela letra 𝐺, chamado de rede recíproca. As redes no espaço real e as redes recíprocas podem ser relacionadas por seus vetores primitivos �̂�1, �̂�2, �̂�3, �̂�1, �̂�2 e �̂�3. Obtendo então: R = α. �̂�1 + β. �̂�2 + γ. �̂�3 e G = α. �̂�1 + β. �̂�2 + γ. �̂�3, sendo que todos vetores da rede recíproca devem satisfazer a seguinte relação: 𝐺. 𝑅 = (𝛼. �̂�1 + 𝛽. �̂�2, +𝛾. �̂�3). (𝛼. �̂�1 + 𝛽. �̂�2, +𝛾. �̂�3) = 2𝜋𝑁 (22) Cada vetor primitivo da rede recíproca pode ser dado por: �̂�1 = 2𝜋 �̂�2 × �̂�3 �̂�1. (�̂�2 × �̂�3) (23) �̂�2 = 2𝜋 �̂�3 × �̂�1 �̂�1. (�̂�2 × �̂�3) (24) �̂�3 = 2𝜋 �̂�1 × �̂�2 �̂�1. (�̂�2 × �̂�3) (25) Portanto, associando os modos magnéticos dos cristais fotônicos ao teorema de Bloch, temos o estado de Bloch: �⃗⃗� 𝑘(𝑟 ) = 𝑒𝑖�⃗� 𝑟 𝑢𝑘(𝑟 ) (26) Onde, �⃗� é um vetor pertencente a G e u(𝑟 ) é a função periódica em um deter- minado sentido de periodicidade e u(𝑟 ) = u(𝑟 + R). Um importante ponto do estado de Bloch é que um dado vetor de onda �⃗� que difere de outro vetor de onda (�⃗� + mb) são idênticos para quaisquer valores de 𝑚 se b = 2π/Λ. Portanto, é necessário somente considerar valores de �⃗� dentro da variação de [-π/Λ à π/Λ]. Essa região é também chamada de Zona de Brillouin. Essa relação 40 permite que uma estrutura periódica seja analisada por meio de uma única célula uni- tária de vetor �⃗� . A Figura 7 ilustra os vetores da rede de Bravais, recíproca e a primeira Zona de Brillouin de um cristal fotônico. Figura 7 – Estrutura de uma cristal fotônico em duas dimensões com arranjo quadrado (a) Rede de Bravais com seus vetores primitivos; (b) Os correspondentes pontos da rede recíproca e seus vetores primitivos; (c) Construção da primeira zona de Brillouin. Fonte: Adaptado de JOANNOPOULOS et al. (1995. p. 237). Portanto, substituindo o teorema de Bloch na Equação Mestre, obtém-se todas as informações referentes ao campo magnético em função do vetor de onda �⃗� . Θ̂𝐻𝑘 = ( 𝜔(�⃗� ) 𝑐 ) 2 𝐻𝑘 (27) ∇ × 1 𝜀(𝑟 ) ∇ × 𝑒𝑖�⃗� .𝑟 𝑢𝑘(𝑟 ) = ( 𝜔(�⃗� ) 𝑐 ) 2 𝑒𝑖�⃗� .𝑟 𝑢𝑘(𝑟 ) (28) (𝑖�⃗� + ∇) × 1 𝜀(𝑟 ) (𝑖�⃗� + ∇) × 𝑢𝑘(𝑟)⃗⃗ ⃗ = ( 𝜔(�⃗� ) 𝑐 ) 2 𝑢𝑘(𝑟 ) (29) Θ̂𝑘𝑢𝑘(𝑟 ) = ( 𝜔(�⃗� ) 𝑐 ) 2 𝑢𝑘(𝑟 ) (30) Onde Θ̂𝑘 é um operador Hermitiano dependente de �⃗� e uk(𝑟 ) = uk(𝑟 + R) = dado por: Θ̂𝑘 = (𝑖�⃗� + ∇) × 1 𝜀(𝑟 ) (𝑖�⃗� + ∇) × (31) 41 Os modos são dados pelos autovalores (ω(�⃗� )/c)2, onde devem atender a condição de transversalidade (i�⃗� + ∇) x uk(𝑟 ) = 0. Para calcular estes modos torna-se necessário a utilização de técnicas computacionais. Essas técnicas não serão abordadas neste estudo, sendo que para as simulações foi utilizado um software de cálculo eletromagnético baseado no método dos elementos finitos (FEM) no domínio da frequência disponibilizado pelos autores (MARRONE; RODRÍGUEZ-ESQUERRE; HERNÁNDEZ-FIGUEROA; 2002) (RODRÍGUEZ-ESQUERRE; HERNÁNDEZ- FIGUEROA; KOSHIBA; 2006). Contudo, há menos de uma década, técnicas utilizando redes neurais artificiais foram propostas para estimar os diagramas de dispersão e a PBG, reduzindo assim o esforço computacional (MALHEIROS- SILVEIRA; HERNANDEZ-FIGUEROA, 2012) (FERREIRA; MALHEIROS-SILVEIRA; HERNÁNDEZ-FIGUEROA, 2018). 3.3 Banda Proibida Fotônica (Photonic Band Gaps - PBG) Os modos eletromagnéticos podem ser classificamos como: transversal elétrico (TE) e transversal magnético (TM). É chamado de modo TE quando o campo elétrico é perpendicular a direção de propagação, e de modo TM quando o campo magnético é perpendicular a direção de propagação. Desta forma, a PBG pode ocorrer tanto para os modos TE quanto TM ou até mesmo para ambos os modos, chamada de banda proibida fotônica absoluta. Cada modo pode estar associado a uma faixa de frequências que atenda a condição de transversalidade, sendo esses modos apresentados de modo discreto e classificados por um índice M, conforme ilustra na Figura 8, onde são apresentados três diagramas de banda fotônica de diferentes combinações de matérias dielétrico em uma estrutura unidimensional. A Figura 8 (a) apresenta o diagrama de banda de uma estrutura homogênea composta de arsenieto de gálio (GaAs), não apresentando nenhuma faixa de frequên- cia proibida. A Figura 8 (b) apresenta uma estrutura composta por arsenieto de gálio e arsenieto de gálio e alumínio (GaAlAs), obtendo uma relação entre as constantes dielétrica maior do que 1(Δɛ > 1). Dessa forma, possui uma faixa de frequência não permitida à propagação de luz, a PBG. Aumentando essa relação entre as constantes dielétricas para Δɛ ≫ 1, a Figura 8 (c) apresenta uma região de PBG ainda maior 42 quando relacionada às demais estruturas. Sendo assim, a PBG ocorre quando ɛ1/ɛ2 ≠ 1. Figura 8 – Estrutura de banda fotônica de cristal fotônico unidimensional com a = 0,5Λ com relação entre as constante dielétrica de: (a) Δɛ = 1; (b) Δɛ = 13/12; (c) Δɛ = 13/1. Fonte: Adaptado de JOANNOPOULOS et al. (1995. p. 46). Conforme foi apresentado, as equações de Maxwell são escaláveis, podendo ser relacionadas a um fator s. Com isso, permite com que as frequências e os vetores de onda sejam apresentados com unidades adimensionais ωΛ/2πc e kΛ/2π, respec- tivamente. Nas quais a relação da estrutura com o comprimento de onda no espaço livre é dada por: λ0 = 2πc/ω. Além disso, a relação da PBG é usualmente denotada mediante à razão chamada de gap-midgap, que condiz com a relação Δω/ωm, sendo Δω a menor espessura entre as duas frequências e ωm a média aritmética entre estas frequências. 3.4 Aplicação de Estruturas Cristalinas Fotônicas Bidimensionais Dentre as diversas aplicações dos cristais fotônicos por meio dos três principais regimes de operação mencionados anteriormente, este trabalho apresentará uma breve revisão das aplicações em regime de reflexão, ocasionada pela PBG. Conforme já mencionado, cristais fotônicos bidimensionais têm a periodicidade em duas dimensões. Um exemplo de tal cristal pode ser representado na Figura 9 (a), composto por colunas cilíndricas de Alumina envolto de ar, em um arranjo quadrado, 43 onde a estrutura é homogênea na direção �̂� e periódico nas direções (�̂�, �̂�), com perí- odo Λ. Na Figura 9 (a) é apresentada uma estrutura cristalina bidimensional composta por colunas cilíndricas de Alumina envolto de ar, onde também demonstra em visão bidimensional (�̂�, �̂�) a relação periódica (Λ) entre as colunas. A estrutura de bandas deste cristal é apresentada na Figura 9 (b), apresentando uma PBG entre os modos TM 1 e 2 (TM12). Figura 9 – (a) Estrutura cristalina com arranjo quadrado de materiais dielétricos em formato de cilindro com raio r e periodicidade Λ, (b) Estrutura de banda fotônica para um arranjo quadrado com colunas de Alumina com n ≈ 2,98, envolto de ar. Fonte: A adaptado de JOANNOPOULOS et al. (1995, p. 67-68). Com base na estrutura apresentada na Figura 9 (a) é possível, por meio de defeitos na periodicidade da estrutura, permitir que se crie uma cavidade ressonante ou um guia de onda com propagação de um ou mais modos, pertencentes a faixa da PBG. 3.4.1 GUIAS DE ONDAS Guias de ondas são utilizados para propagar ondas de um ponto a outro, à uma determinada frequência, geralmente confinando a luz por meio da diferença de índices de refração entre a parede e o guia, propriamente dito, dada pela Lei de Snell. Outro método de propagar da luz pode ocorrer confinando-a numa estrutura fotônica por meio da sua PBG. Onde, ocasionando um defeito ou removendo ou alte- rando o diâmetro de uma ou mais linhas colunas da estrutura apresentada na Figura 9 (a) permite com que uma ou mais frequências pertencentes a sua PBG, se propague no seu interior. 44 A Figura 10 apresenta dois exemplos de estrutura cristalina composta por lacu- nas de ar numa superfície de silício que atuam como guias de ondas, onde na Figura 10 (a) é ocasionado um defeito em uma linha de lacunas destacada em vermelho, preenchendo-a de silício, ocasionando um guia de onda monomodo conforme ilustra o diagrama de banda, e na Figura 10 (b) o defeito é ocasionado em duas linhas des- tacadas em vermelho do arranjo periódico, onde permite com que dois modos perten- centes a PBG se propaguem entre a estrutura, conforme ilustra o seu diagrama de banda. Figura 10 - Guia de onda reto em uma estrutura cristalina fotônica bidimensional. (a) diagrama de banda para os vetores de onda na mesma direção da linha com defeito e (b) a propagação da luz com a polarização do campo �⃗� 𝑧. Fonte: Adaptado de BUSCH et al. (2006). Além do guia de onda reto, a propriedade de confinar a luz por meio da PBG permite com que o guia de onda seja projetado com curva de 90º. A Figura 11 ilustra a imagem microscopia de um exemplo de tal guia de onda otimizado e projetado num cristal fotônico bidimensional com arranjo triangular por (CHANG-ZHU; YA-ZHAO; ZHI-YUAN, 2010). Figura 11 – Imagem microscopia de um guia de onda com curva de 90º em uma estrutura cristalina fotônica. Fonte: (CHANG-ZHU; YA-ZHAO; ZHI-YUAN, 2010). 45 3.4.2 FILTRO PASSA FAIXA Um defeito pontual como, por exemplo, a remoção ou alteração do raio de uma única coluna da estrutura apresentada na Figura 9 (a), ocasiona o efeito de cavidade ressonante para uma determinada frequência pertencente a faixa de PBG, podendo ser utilizado para confecção de laser ou filtros. Um exemplo de filtro ocorrer quando unimos dois guias de onda por meio de uma cavidade ressonante, tornando-se possível projetar uma estrutura que propague a luz, de um guia ao outro, somente para a frequência de ressonância. Um exemplo desta estrutura foi projetado por (ZHOU et al., 2014) em um cristal fotônico bidimensional com arranjos triangulares formados por lacunas de ar numa superfície de silício conforme ilustra a Figura 12 (a). Tal estrutura permite a propagação da luz em somente uma pequena faixa de frequência conforme demostra o seu gráfico de transmissão em relação às frequências é apresentado na Figura 12 (b). Figura 12 – Exemplo de filtro passa faixa, (a) estrutura cristalina e suas respectivas relações geométricas (b) gráfico de transmissão em relação as frequências para relações:r1 = 0,29Λ e r2 = 0,18Λ. Fonte: Adaptado de ZHOU et al. (2014). 3.4.3 DIVISOR DE POTÊNCIA Outro possível dispositivo explorando a PBG dos cristais fotônicos bidimensionais é o divisor de potência, formado por um guia de onda com uma entrada e duas ou mais saídas. A Figura 13 ilustra um exemplo desta aplicação em uma estrutura cristalina de colunas cilíndricas de dielétrico envolto de ar, contendo um guia de onda com uma entrada e duas saídas. 46 Figura 13 - Divisor de potência em cristal fotônico bidimensional. (a) Esquemático da estrutura formado por linhas de defeitos reduzindo o diâmetro das colunas de dielétrico; (b) Propagação do campo �⃗� 𝑧 dentre a estrutura. Fonte: Adaptado de LIU (2012). 4 METODOLOGIA Inicialmente, foi realizado um estudo comparativo entre o ABC e o GA por meio de duas funções multimodais previamente conhecidas, com o objetivo de compreender melhor o comportamento de busca do ABC, para posteriormente aplicá- lo à otimização de estruturas cristalinas fotônicas bidimensionais em arranjo triangular. Todas as simulações deste trabalho foram executadas na plataforma Matlab® em um notebook com processador Intel® i7 2,7Ghz, 8 Gb de RAM e Windows 10. 4.1 Comparação entre ABC e GA em Funções Multimodais Continuas. Utilizando funções previamente conhecidas, nas quais são conhecidos, antecipadamente, os seus valores mínimos e máximos, torna-se possível avaliar a eficiência quanto aos resultados alcançados por cada meta-heurística, além de compreender melhor o seu comportamento de busca, por meio de gráficos. As funções utilizadas e os intervalos de seus domínios são apresentados abaixo, sendo que a função objetivo foi estabelecida pela maximização de cada uma dessas funções. Função 1: 𝑓1(𝑥, 𝑦) = 0,5 + 𝑠𝑒𝑛2(√𝑥2 + 𝑦2) − 0,5 (1 + 0,001(𝑥2 + 𝑦2)) 𝑥, 𝑦 ∈ [−10,10] Função 2: 𝑓2(𝑥, 𝑦) = 𝑥. 𝑠𝑒𝑛(4. 𝜋. 𝑥) − 𝑦. 𝑠𝑒𝑛(4. 𝜋. 𝑦 + 𝜋) + 1 𝑥, 𝑦 ∈ [−2,2] 47 Ambos os algoritmos empregados são populacionais e possuem um critério de parada, portanto, os parâmetros denominados critério de parada e população total foram igualmente configurados para as duas meta-heurísticas. O critério de parada foi estabelecido pela quantidade máxima de ciclos, sendo simulado para quatro diferentes valores: 100, 200, 400 e 800 ciclos. Para facilitar a comparação, cada solução do GA (cromossomos) e do ABC (abelhas), será tratada, igualmente neste trabalho, como um indivíduo. Nas simulações, foram estabelecidos quatro diferentes valores para a quantidade de indivíduo em uma população, sendo: 10, 20, 40 e 80 indivíduos. Para cada uma das funções foram realizadas 16 diferentes simulações correspondentes a cada possíveis combinações dentre os valores propostos acima para os paramentos “critério de parada” e “tamanho da população”. Estas variações, entre a quantidade de população e o número máximo de ciclos, foram propostas com a finalidade de compreender as suas implicações no resultado final. Para garantir uma melhor eficácia dos resultados e eliminar resultados gerados ocasionalmente, cada função foi submetida a um número de 10 simulações diferentes para cada conjunto de configurações, sendo apresentado como o ótimo de cada função o resultado da média aritmética dos valores encontrados. Outros parâmetros de controle do GA e ABC foram igualmente configurados para todas as 16 simulações. A Tabela 1 relaciona cada um deles com as atribuições adotadas. Tabela 1 – Parâmetros de controle do GA e ABC para otimização das funções bidimensionais. Algoritmo Parâmetros Atribuição GA Cromossomos: Cada variável é presentada por um número real de 16 bits Tipo de Seleção: Roleta Taxa de Mutação: 1% alterando 2 cromossomos escolhido aleatoriamente a cada geração. Taxa de Crossover: 95% com 1 posição de corte aleatória. ABC Abelhas Empregadas: 50% da População total Abelhas Espectadoras: 50% da População total Parâmetro limite 4 ciclos No ABC, cada abelha é representada por dois números reais equivalente as duas variáveis das Funções (x,y). Para cada simulação foi gerada, aleatoriamente, 48 uma única população inicial para ambos os algoritmos. A equação (31) foi utilizada para gerar cada indivíduo da população inicial. 𝑝𝑜𝑝𝑖(𝑥, 𝑦) = (𝑟𝑎𝑛𝑑(0,1) . (𝑚𝑎𝑥 − 𝑚𝑖𝑛)) + 𝑚𝑖𝑛 (31) Onde: pop é um vetor contendo os valores de cada variável x e y do indivíduo de índice i. rand(0,1) é um número real aleatório entre [0,1]. max e min são os valores máximo e mínimo do domínio de cada variável da função. 4.1.1 SIMULAÇÃO COMPARATIVA DA FUNÇÃO 1 A Função 1 tem o comportamento apresentado na Figura 14, contendo não somente um ponto de máximo, mas diferentes regiões em formato de anéis. O anel central é a região que apresenta o maior valor da função, f1(x,y) ≈ 0,99877, sendo então a região pertencente aos ótimos globais da função. Os anéis subsequentes apresentam ótimos com diferenças menores do que um centésimo em relação ao ótimo global. A diversidade de regiões contendo ótimos locais com pequena diferença auxilia na comparação da eficiência do método de busca entre cada meta-heurística. Figura 14 – Gráfico de superfície Função 1. A média dos ótimos encontrada por cada algoritmo é apresentada na Figura 15, estruturada em quatro gráficos diferentes para cada quantidade de população proposta. Cada gráfico em barra compara os resultados obtidos pelos algoritmos para diferentes configurações de parada. 49 O ABC apresentou resultados superiores em todas as configurações propostas quando comparado ao GA, tendo o pior resultado para a configuração de 10 indivíduos e 100 ciclos, obtendo um resultado de 0,2% menor do que o máximo global da função. Além disso, o ABC apresentou, em todas as configurações cujos número de ciclos foram maiores ou iguais a 200, o valor do ótimo global da função, f1(x,y) ≈ 0,99877. Analisando mais a fundo, a Figura 16 apresenta a posição final de cada um dos indivíduos de ambos os algoritmos para configuração de 80 indivíduos e 800 ciclos. Sendo representadas em vermelho todas as abelhas empregadas do ABC e em amarelo todos os cromossomos do GA. É possível notar no gráfico da Figura 16 uma maior densidade de indivíduos do ABC na região do anel central, onde contém os ótimos globais da função. Isso mostra, para as configurações adotadas para cada algoritmo, que o ABC conseguiu obter uma diversidade de ótimos locais. Já o GA tem a característica de convergir a maioria dos indivíduos para o mesmo ponto, limitando-se no espaço de busca. Um dos parâmetros que pode modificar essa característica do GA é a taxa de mutação, aumentando assim sua diversidade de busca. Figura 15 – Resultados das simulações referente a Função 1 com todos os parâmetros de parada propostos para (a) 10 indivíduos, (b) 20 indivíduos, (c) 40 indivíduos e (d) 80 indivíduos. 50 Figura 16 – Gráfico tridimensional da função 1 com as posições finais de cada indivíduo para 80 indivíduos e 800 ciclos. Para analisar melhor o espaço de busca dos algoritmos, foram gerados dois gráficos contendo a evolução de cada indivíduo no decorrer de toda a simulação, sendo ilustrado na Figura 17 uma relação de ótimo encontrado a cada ciclo. É possível notar que ambos os algoritmos obtiveram uma convergência rápida para o ótimo encontrado, porém o ABC aparenta uma maior exploração do espaço de busca em relação ao GA. 4.1.2 ANÁLISE COMPARATIVA DA FUNÇÃO 2 A Função 2 apresenta ótimos locais ao longo de todo o seu intervalo conforme ilustra a Figura 18, contendo quatro distintas regiões com os ótimos globais, f2(x,y) ≈ 4,2539. Essas regiões se encontram próximas as quatro extremidades do gráfico da Figura 18. Figura 17 – Gráfico de evolução para Função 1 com configuração de 80 indivíduos e 800 ciclos dos algoritmos. (a) Algoritmo ABC e (b) Algoritmo GA. (a) (b) 51 Todas análises realizadas para a Função 1 foram igualmente propostas para a Função 2. Os resultados da média entre as simulações para cada conjunto de configuração estabelecido à Função 2, são apresentados na Figura 19. Figura 18 – Gráfico de superfície Função 2. Para Função 2, o GA apresentou um melhor desempenho para as configurações do parâmetro de parada com 100 ciclos, demostrando que a dinâmica de busca do GA permite encontrar um ótimo local melhor do que o ABC com o número Figura 19 – Resultados das simulações para função 2 para cada configuração, dividindo em quatro gráficos para diferentes quantidades de população (a) 10 indivíduos, (b) 20 indivíduos, (c) 40 indivíduos e (d) 80 indivíduos. 52 de indivíduos menor e poucos ciclos, porém, o resultado alcançado está relativamente longe do ótimo global da função. O ABC alcançou o resultado igual ao ótimo global em grande parte das configurações, já o GA não alcançou em nenhuma das configurações, obtendo seu melhor resultado com a configuração de 40 indivíduos e 400 ciclos, equivalente a 0,6% inferior em relação ao valor máximo da função. A Figura 20 mostra o gráfico da Função com a posição final de cada indivíduo para a configuração de 80 indivíduos e 800 ciclos. É possível notar a diversidade no espaço de busca do ABC em relação ao GA, onde os indivíduos do ABC estão divididos em três regiões distintas, as quais contém os ótimos globais da função. Já no GA ocorre a convergência de todos os indivíduos para somente um ponto máximo da função. Observando também o gráfico de evolução ilustrado na Figura 21, o ABC apresentou a convergência ao melhor resultado com uma quantidade menor de ciclos em comparação ao GA. Além disso, o ABC obteve uma maior diversidade dos valores da função objetivo, mostrando uma maior área de busca no domínio da função. Figura 20 – Gráfico tridimensional da função 2 com as posições finais de cada indivíduo para 80 indivíduos e 800 ciclos. 53 Figura 21 – Gráfico de evolução para Função 1 com configuração de 80 indivíduos e 800 ciclos dos algoritmos. (a) Algoritmo ABC e (b) Algoritmo GA. Para as mesmas simulações da análise anterior, com 80 indivíduos e 800 ciclos, foi proposta uma análise do recurso computacional utilizado, mensurando o tempo utilizado para concluir toda a simulação de cada algoritmo, para ambas as Funções, conforme apresentado na Figura 22. O ABC operou com um tempo inferior ao GA de aproximadamente 41,2% para a Função 1 e 36,31% para a Função 2. Um dos fatores que impactam nessa diferença é que a interação do ABC com a função objetivo ocorre, em grande parte, somente pela quantidade de abelhas espectadoras mais uma abelha exploradora a cada ciclo, já o GA é pela quantidade de população total a cada ciclo. Portanto, através dessas análises, foi possível compreender melhor o comportamento dos algoritmos para uma implementação de otimização real. O ABC apresenta melhores resultados em comparação ao GA, tanto nos ótimos encontrados quanto no tempo computacional utilizado, para ambas funções, além de demostrar (a) (b) Figura 22 – Média do tempo computacional total utilizado pelo GA e o ABC, para configuração de 80 indivíduos e 800 ciclos, para: (a) Função 1 e (b) Função 2. 54 que utiliza um menor recurso computacional pelo fato de interagir com a função objetivo um número igual a quantidade de abelhas espectadoras mais uma abelha exploradora. 4.2 Otimização de Estrutura Cristalina Fotônica Bidimensional Na Figura 9 (a) foi apresentado um exemplo de design de um cristal fotônico bidimensional composto por colunas de alumina envoltas de ar. A Figura 23 apresenta, em visão bidimensional, quatro exemplos de design, nos quais o material dielétrico com maior índice de refração é representado pela cor preta e com menor índice na cor branca. A configuração da Figura 23 (b) corresponde ao mesmo arranjo apresentado na Figura 9 (a). Figura 23 – Visão bidimensional de estruturas cristalinas fotônicas sendo formada por: (a) Lacunas cilíndricas de ar em arranjo quadrado; (b) Cilindros dielétricos em arranjo quadrado; (c) Lacunas cilíndricas de ar em arranjo triangular; (d) Cilindros dielétricos em arranjo triangular. A PBG, para os modos TM, tende a ser maximizada em estruturas de colunas dielétricas envoltas em ar. Já para os modos TE, em estruturas com colunas de ar em uma placa dielétrica. Além disso, estruturas com arranjo triangular apresentam uma maior faixa de PBG quando comparadas àquelas com arranjos quadrado (JOANNOPOULOS et al., 1995). 4.3 Adaptação do ABC para Otimização do Cristal Fotônico Neste trabalho foram utilizadas estruturas em arranjo triangular compostas de Si com índice de refração ne = 3,477 e ar com n0 = 1, aplicando o ABC a maximizar a PBG entre os modos TE12, TE34, TM12 e TM34 em quatro diferentes estruturas. Os 55 resultados obtidos foram comparados aos obtidos por estruturas otimizadas pelo GA em (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007). O problema foi tratado como a busca inversa de um novo design de uma célula unitária, utilizando, para o cálculo da estrutura de banda, um software de cálculo eletromagnético baseado no método dos elementos finitos (FEM) no domínio da frequência (MARRONE; RODRÍGUEZ-ESQUERRE; HERNÁNDEZ-FIGUEROA; 2002). A função objetivo é dada pela equação (32). PBG(%) = 𝐸𝑐𝑖𝑚𝑎 − 𝐸𝑏𝑎𝑖𝑥𝑜 𝐸𝑚𝑒𝑑𝑖𝑎 . 100% (32) Onde Ecima é a menor frequência do modo de ordem mais alta, Ebaixo é a maior frequência do modo de ordem mais baixa e Emedia a média dos dois valores. A célula unitária foi discretizada em um conjunto de 100 elementos enumerados, formando um arranjo triangular de 10x10. A Figura 24 ilustra um exemplo desta estrutura com 9 elementos representada em: (a) um cristal fotônico com 16 células unitárias; (b) uma célula unitária e (c) a representação computacional de uma abelha. Cada elemento da célula unitária é representando, computacionalmente, por um número binário e representado pelo índice j de uma abelha, sendo que a otimização da estrutura constitui em uma busca discreta, cuja a configuração torna possível 2100 soluções. O ABC foi apresentado por Karaboga para otimização de funções contínuas (KARABOGA; BASTURK, 2008) porém, para otimizar estrutura fotônica em um domínio discreto, fez-se necessário algumas alterações da versão clássica, sendo elas: Figura 24 – (a) Estrutura com arranjo triangular cristalina composta por 16 células unitárias; (b) Célula unitária representado por 9 elementos e (c) Representação de uma abelha em modo binário para célula unitária ao lado. 56 • Cada abelha no ABC foi representada por um vetor de 100 posições (j = 100), sendo que cada posição pode assumir somente valores binários, 0 ou 1, indicando a ausência ou presença de silício, respectivamente. • Para aumentar a diversidade, somente posições escolhidas aleatoriamente são sujeitas à “dança”, as demais são atribuídas os valores correspondentes a abelha empregada escolhida para a “dança”. • A “dança das abelhas” representada pela equação (2), foi adaptada para a lógica descrita no Algoritmo 2, garantindo que cada posição do vetor contenha apenas valores 0 ou 1. Onde x e y são os valores correspondentes as abelhas empregada e espectadora, respectivamente, de índices i, representando a i-enésima abelha, j é a posição de seu vetor, k o índice aleatório de uma abelha empregada diferente de i e a variável P correspondente a um número binário aleatório utilizado para determinar qual posição será sujeita a “dança”, ou não, sendo P = 1 o elemento de posição j é sujeita a “dança” e quando P = 0 é atribuído o valor da abelha empregada selecionada, conforme mencionado no item acima. Algoritmo 2 – Adaptação do ABC para o domínio discreto 1- 𝑃𝑗 = 𝑎𝑙𝑒𝑎𝑡𝑜𝑟𝑖𝑜(0 𝑜𝑢 1) 2- SE 𝑃𝑗 = 1 FAZER: 3- 𝑦𝑖𝑗 = 𝑥𝑖𝑗 + 𝜃 ∗ ‖𝑥𝑖𝑗 − 𝑥𝑘𝑗‖; 4- SE 𝑦𝑖𝑗 > 0 FAZER: 5- 𝑦𝑖𝑗 = 1; 6- SENÃO: 7- 𝑦𝑖𝑗 = 0; 8- FIM; 9- SENÃO: 10- 𝑦𝑖𝑗 = 𝑥𝑖𝑗 11- FIM; Os parâmetros utilizados para otimização dos cristais são apresentados na Ta- bela 2, respeitando os mesmos parâmetros do ABC clássico. A população inicial foi gerada de modo aleatório, para todas as simulações. Para cada modo foram realiza- das uma simulação diferente, tendo como resultado quatro diferentes estruturas, uma para cada modo. 57 Tabela 2 – Parâmetros de controle ABC para otimização do cristal fotônico. Parâmetros Atribuição População total: 40 abelhas Abelhas empregadas: 50% da População total Abelhas espectadoras: 50% da População total Parâmetro limite 10 ciclos Máximo número de ciclos: 2000 A quantidade total de abelhas do ABC foi selecionada com base nos experi- mentos computacionais realizados entre as duas funções multimodais, nas quais no ABC a interação com a função objetivo para busca por novas fontes de alimentos ocorre de acordo com a quantidade de abelhas espectadoras mais uma abelha explo- radora a cada ciclo. Sendo assim, para comparar com os resultados alcançados no GA, foi pro- posto, neste trabalho, uma quantidade de abelhas espectadoras semelhante a quan- tidade de cromossomos utilizado pelo GA em MALHEIROS-SILVEIRA; RODRÍGUEZ- ESQUERRE, 2007. Portanto, a população total do ABC, constituída por 50% de abe- lhas espectadoras e 50% abelhas empregadas torna-se igual ao dobro da população utilizada pelo GA, proporcionando uma comparação mais justa entre ambos os algo- ritmos quando comparado pela quantidade de interações com a função objetivo. Con- tudo, após todas simulações, o ABC foi submetido ao processo de otimização de 4 novas estruturas, uma para cada modo proposto, com o número total de abelhas igual a 20. 5 RESULTADOS Em (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007), os autores realizaram três diferentes testes, modificando alguns parâmetros do GA, sendo eles: (1) Taxa de mutação, taxa de cruzamento e método de seleção; (2) Tamanho da população; (3) Os métodos de inicialização da população. A Tabela 3 apresenta os parâmetros que apresentaram os melhores resultados, sendo a taxa de cruzamento e o ponto de cruzamento determinados aleatoriamente, assim como a taxa de mutação (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007). 58 Tabela 3 – Parâmetros adotados para otimização com GA. Parâmetros Atribuição População total: 20 indivíduos Máximo número de ciclos: 2000 Tipo de Seleção: Roleta Tamanho dos cromossomos: 100 bits Os autores submeteram o algoritmo a fim de otimizar uma estrutura para o modo TE12 e outra para TM34, obtendo as células unitárias apresentadas na Figura 25. As PBGs encontradas pelo GA foram de 47,85% para o modo TE12 e 34,31% para o modo TM34. Entretanto, o índice de refração utilizado para o Silício foi de ne = 3,4. Já para este trabalho foi utilizado um valor mais atualizado de ne = 3,477 (SALZBERG; VILLA, 1957). Portanto, foram recalculados as PBG’s para cada modo com o novo índice, obtendo os valores de PBG de 49,11% para o modo TE12 e 34,78% para o modo TM34. Fonte: (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007). Cabe destacar que a estrutura e o valor da PBG apresentados por (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007), referente ao modo TM12, estão divergentes, sendo o modo correto o TE12, tendo ocorrido um erro de digitação. Para testar as alterações do ABC foi realizada uma simulação com 500 ciclos e 20 abelhas, maximizando a PBG entre os modos TE12. O resultado encontrado foi de uma PBG de 46,95%, conforme ilustra o diagrama de bandas na Figura 26. Figura 25 – Cristal fotônicos otimizada pelo GA (a) para o modo TE12 e (b) para o Modo TM34. 59 Figura 26 – Diagrama de banda de estrutura teste para o modo TE12. O resultado encontrado foi próximo ao da estrutura otimizada pelo GA de 49,11%, sendo que o ABC utilizou aproximadamente menos de um terço dos recursos computacionais, convergindo para o melhor resultado próximo ao ciclo 400. O gráfico de evolução do ABC para otimização e a célula unitária da estrutura encontrada são apresentados na Figura 27, sendo que os pontos em preto são os valores da variável Gbest em cada ciclos e em vermelho a PBG de cada abelha empregada. Figura 27 – (a) Gráfico de evolução do ABC para o modo TE12 com 20 abelhas e 500 ciclos e (b) a célula unitária encontrada com PBG de 45,59%. 5.1 Otimização com População Total de Abelhas Igual a 40 Validadas as alterações do algoritmo, o ABC foi submetido à otimização para ambos os modos otimizados pelo GA, modos TE12 e TM34. O ABC apresentou uma PBG(%) de 49,11% para o modo TE12 e 43,31% para o TM34, conforme ilustram os diagramas de banda nas Figura 28 e Figura 29, respectivamente. 60 Figura 28 – Diagrama de banda de estrutura para o modo TE12. Figura 29 – Diagrama de banda de estrutura para o modo TM34. A Figura 30 apresenta os cristais fotônicos e as células unitárias encontradas. Conforme mencionado anteriormente, os modos TM têm a PBG maximizada para estruturas com colunas dielétricas envoltas em ar e os modos TE em estruturas com colunas de ar, sendo, portanto, os resultados encontrados coerentes com a literatura (JOANNOPOULOS et al., 1995). Figura 30 – Cristal fotônicas e célula unitária otimizada pelo ABC: (a) Modo TE12; (b) Modo TM34. O valor da PBG da estrutura otimizada pelo ABC para o modo TE12 foi igual ao encontrado pelo GA, já para o modo TM34, houve um acréscimo de 24,5% na PBG quando comparada a estrutura otimizada pelo GA. Analisando a curva de evolução do GA para o modo TE12 apresentada no gráfico superior da Figura 31, nota-se que a convergência para o melhor resultado encontrado ocorre entorno da geração 1300. Já a curva de evolução do ABC, ilustrada 61 no gráfico inferior na Figura 31, que apresenta em preto o Gbest e em vermelho a quantidade de néctar em cada fonte de alimento das abelhas empregadas, demonstra uma convergência ao melhor resultado encontrado próximo a 200 ciclos, além de encontrar uma maior diversidade de PBG’s no decorrer de 2000 ciclos. Figura 31 – Gráfico de evolução do GA e ABC para o modo TE12. Fonte: (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007). Analogamente, para o modo TM34, a curva de evolução do GA ilustrada no gráfico superior na Figura 32, apresenta uma convergência ao melhor resultado encontrado próximo a 1300 a 1400 gerações. Já o gráfico inferior na Figura 32 apresenta a curva de evolução do ABC para o modo TM34 e mostra uma convergência mais rápida em relação ao GA para o melhor resultado encontrado, próximo ao ciclo 400. 62 Figura 32 – Gráfico de evolução do GA e ABC para o modo TM34. Fonte: (MALHEIROS-SILVEIRA; RODRÍGUEZ-ESQUERRE, 2007). A Tabela 4 apresenta um comparativo com os respectivos ganhos entre os resultados obtidos pelo ABC para cada modo. Tabela 4 - Comparação de resultados entre GA e ABC. Modo PBG (%) Relação de Ganho ABC GA TE12 49,11% 49,11% +0,0% TM34 43,31% 34,78% +24,5% Além dos modos otimizados pelo GA, foram também realizadas a otimização de estruturas para os modos TM12 e TE34. A Figura 33 apresenta o diagrama de banda para o modo TM12, alcançando uma PBG de 47,55%. 63 Figura 33 – Diagrama de banda de estrutura para o modo TM12. Para o modo TE34, a PBG encontrada foi de 10,37%, conforme ilustra a Figura 34, 10,37%. Figura 34 – Diagrama de banda de estrutura para o modo TE34. Os cristais fotônicos e suas respectivas células unitárias, referentes aos ótimos encontrados nos casos representados nas Figuras 33 e 34, são apresentados na Figura 35. Figura 35 – Cristais fotônicos e célula unitária otimizada pelo ABC: (a) Modo TM12; (b) Modo TE34. O gráfico de evolução para otimização do modo TM12 é apresentado na Figura 36, apresentando uma convergência ao melhor resultado próximo a 400 ciclos. 64 Figura 36 – Gráfico de evolução do ABC para o modo TM12. Para todos os modos foram considerados uma condição mínima de espaço entre os diferentes dielétricos, tanto na horizontal, vertical e diagonal de um elemento, equivalendo a um décimo do comprimento do período Λ, garantindo assim que as estruturas sejam compatíveis com a fabricação em dispositivos fotônicos integrados. Portanto, para o modo TE34, foram descartados resultados maiores, entorno de 17%, que não atendiam esse critério. A Figura 37 apresenta o gráfico de evolução para o modo TE34, desconsiderando todos os resultados descartados, por não atenderem o espaçamento mínimo entre os dielétricos. Figura 37 – Gráfico de evolução do ABC para o modo TE34. 5.2 Otimização com População Total de Abelhas Igual a 20 Submetido à otimizar quatro novas estruturas com a quantiade total de abelhas igual a 20, respeitando a mesma proporção de 50% abelhas empregadas e 50% abelhas espectadoras e o parametro limite igual a 10, o ABC alcançou os resultados apresentados na Tabela 5, para cada modo proposto. 65 Tabela 5 – Resultados obtidos na otimização com ABC com quantidade de abelhas totais igual 20. Modo PBG TM12 47,06% TM34 40,99% TE12 49,20% TE34 9,99% Os respectivos cristais fotônicos formados pelas células unitárias obtidos na simulação com o número de abelhas reduzido (20) estão ilustrados na Figura 38. Figura 38 – Cristais fotônicos compondo as células unitárias encontrada na otimização com quantidade de abelhas totais igual a 20 para os seguintes modos: (a) TE12; (b) TE34; (c) TM12 e (d) TM34. A Figura 39 apresenta os gráficos de evolução de cada um dos modos otimizados, apresentando somente a evolução até 1000 ciclos. Contudo, ambas as simulações decorreram 2000 ciclos, porém o ABC se estabilizou antes, conforme ocorreu para as otimizações com a quantidade total de abelhas igual a 40. 66 Figura 39 – Gráficos de evolução para os quatros modos propostos com a quantidade de abelhas totais igual a 20. Fazendo uma análise comparativa entre a evolução do ABC, tanto com 20 quanto com 40 abelhas em relação ao GA, podemos notar que, o ABC obteve uma convergência ao máximo valor encontrado de PBG mais rápida para ambas as configurações propostas, conforme demonstram as Figuras 40 e 41, para os modos TE12 e TM34, respectivamente. Figura 40 – Gráfico comparativo entre a evolução do ABC e do AG para maximização da PBG do modo TE12. 67 Figura 41 – Gráfico comparativo entre a evolução do ABC e do AG para maximização da PBG do modo TM34. Para as simulações com 20 abelhas para os modos TE12 e TM34, foi mensurado a quantidade de vezes que o ABC faz chamada ao FEM, considerando que o maior esforço computacional do algoritmo consiste no cálculo da estrutura de banda de cada uma das estruturas encontrada, onde obteve para o moto TE12 20.545 chamadas ao FEM e para o modo TM34 20.748 chamadas ao FEM. Como no GA proposto por MALHEIROS-SILVEIRA; RODRÍGUEZ- ESQUERRE, 2007 foi estabelecido uma população estática a cada geração, é possível mensurar a quantidade mínima que o algoritmo fez chamada à função fitness, por meio do produto da quantidade de indivíduos (20) pelo número de gerações (2000). Sendo assim, o GA fez chamada ao FEM no mínimo 40.000 vezes. Contudo, como ambas as otimizações utilizam o mesmo software para calcular a estrutura de banda de cada celula unitária, podemos notar que o ABC com 20 abelhas utilizou 48% a menos do recurso computacional quando comparado ao GA, obtendo resultados de PBG superiores. Por fim, os resultados obtidos pelo ABC, em todas as simulações realizadas, são apresentados na Tabela 6. Ocorrendo que para o modo TE12, o ABC alcançou um resultado superior com a quantidade de 20 abelhas em comparação com a configuração de 40 abelhas, isso ocorre pelo fato de o ABC ser um algoritmo probabilístico, não garantindo o máximo global a cada simulação, mas sim a aproximação de tal. 68 Tabela 6 – Resultados das PBGs obtidas por estruturas otimizadas pelo ABC. Modos PBG 40 Abelhas 20 Abelhas TM12 47,55% 47,06% TM34 43,31% 40,99% TE12 49,11% 49,20% TE34 10,37% 9,99% 6 CONSIDERAÇÕES FINAIS Nos últimos anos, o ABC tem se demonstrado eficiente na otimização de alguns problemas de engenharia. Neste trabalho foi proposto, pela primeira vez, a utilização do ABC para otimização inversa de cristais fotônicos para maximizar a PBG para os modos TE e TM. Foram necessárias algumas alterações para resolução de problemas em um domínio discreto, haja visto que o ABC clássico foi proposto a otimizar somente função contínua em todo o seu domínio. Porém, as alterações foram cuidadosamente realizadas, buscando o mínimo de alteração possível, garantindo assim a simplicidade do algoritmo apresentado por Karaboga. Inicialmente neste trabalho foi proposto a comparação entre o ABC e o GA por meio de duas funções multimodais, obtendo uma conclusão satisfatória ao ABC, quanto aos resultados e quanto ao compreendimento da dinâmica de busca do ABC. Sendo observado que o esforço computacional do ABC, quando comparado com a interação da função objetiv