## Abstract

We investigate the use of a Genetic Algorithm (GA) to design a set of photonic crystals (PCs) in one and two dimensions. Our flexible design methodology allows us to optimize PC structures for specific objectives. In this paper, we report the results of several such GA-based PC optimizations. We show that the GA performs well even in very complex design spaces, and therefore has great potential as a robust design tool in a range of PC applications.

©2007 Optical Society of America

## 1. Introduction

Photonic crystals (PCs) describe a class of semiconductor structures with a periodic variation of refractive index in 1, 2, or 3 dimensions. As a result, PCs possess a photonic band gap – a range of frequencies in which the propagation of light is forbidden [1, 2]. This unique characteristic of PCs enables them to be used to manipulate light. PCs have already been used for applications such as modifying the spontaneous emission rate of emitters [3, 4], slowing down the group velocity of light [5, 6], and designing highly efficient nanoscale lasers [7].

Given that photonic crystals find applications in a myriad of areas, we proceed to investigate the question: *What is the best possible manufacturable PC design for a given application?* Man-ufacturability refers to the ability to realize the structure with standard tools such as electron beam lithography. Traditionally, photonic crystal designs have been optimized largely by either trial-and-error, iterative searches through a design space, by physical intuition, or some combination of the above methods [8, 9]. However, such methods of design have their limitations, and recent developments in PC design optimization have instead taken on a more systematic and algorithmic nature [10, 11, 12, 13]. In this work, we report the results of a Genetic Algorithm (GA) to optimize the design of a set of one and two-dimensional PC structures. We show that the GA can effectively optimize PC structures for any given design objective, and is thus a highly robust and useful design tool.

## 2. Genetic algorithms

Genetic algorithms (also known as Evolutionary algorithms) are a class of optimization algorithms that apply principles of natural evolution to optimize a given objective [14, 15, 16]. In the genetic optimization of a problem, different solutions to the problem are picked (usually randomly), and a measure of fitness is assigned to each solution. On a given generation of the design, a set of operations, analogous to mutation and reproduction in natural selection, are performed on these solutions to create a new generation of solutions, which should theoretically be “fitter” than their parents. This process is repeated until the algorithm terminates, typically after a pre-defined number of generations, or after a particularly “fit” solution is found, or more generally, when a generation of solutions meets some pre-defined convergence criterion. In general, GA’s are best suited to problems which have a low computational cost for each fitness evaluation.

## 3. Implementation

Genetic algorithms have already been used in PC design to find non-intuitive large-bandgap designs [12, 17] and for designing PC fibers [18]. In this work, we implement a general GA to optimize 1 and 2-dimensional photonic bandgap structures, and show that it is able to robustly optimize these structures for a wide variety of objectives. In the 1-dimensional case, we consider the design of planar photonic crystal cavities by varying the widths of dielectric stacks, accounting for the remaining spatial dimensions via an effective index model; in the 2-dimensional case, we perform the genetic optimization by varying the sizes of circular holes in a triangular lattice. These approaches were chosen because the search space is conveniently well-constrained in these paradigms, and the optimized structures (for the triangular lattice) can be easily fabricated.

In our genetic algorithm, the variables to be optimized (stack widths in 1D, hole radii in 2D) were directly stored in a vector, called a chromosome. Each chromosome therefore compactly represents a dielectric structure to be simulated. Our implementation can be easily modified to optimize over other types of variables as well, such as the positions of the various holes, or the refractive index of the dielectric material.

To propagate a new generation of chromosomes from the current generation, we used the following steps:

**Selection.**

We used fitness-proportionate selection (also known as roulette-wheel selection), to choose parent chromosomes for mating. In this selection scheme, a chromosome is selected out of N chromosomes with a probability *P _{i}* that is proportional to its fitness

*f*, as shown in Eq. (1).

_{i}**Mating.**

After a pair of parent chromosomes *v _{parent}*,1 and

*v*,2 were selected, they were mated to produce a child chromosome

_{parent}*v*by taking a random convex combination of the parent vectors, as in Eq. (2).

_{child}$${\overrightarrow{v}}_{\mathrm{child}}=\lambda {\overrightarrow{v}}_{\mathrm{parent},1}+\left(1-\lambda \right){\overrightarrow{v}}_{\mathrm{parent},2}$$

**Mutation.**

Mutation was used to introduce diversity in the population. We used two types of mutation in our simulations, a random-point crossover and a Gaussian mutation.

**Random-point crossover**: For an original chromosome vector*v*→_{orig}of length*N*, we selected a random index,*k*, from 0 to*N*- 1 as the crossover point, and swapped the two halves of*v*→_{orig}to produce the mutated vector,*v*→_{mut}, as represented in Eq. (3).**Gaussian mutation**: To mutate a chromosome vector by Gaussian mutation, we define each element of*v*→_{mut}to be independent and identically distributed Gaussian Random Variables with mean*v*→_{orig}and standard deviation of σ. This searches the space in the vicinity of the original chromosome vector*v*→_{orig}.$${v}_{i}^{\mathrm{mut}}~N\left({v}_{i}^{\mathrm{orig}},{\sigma}^{2}\right)\phantom{\rule{0.8em}{0ex}}i\in \{1,2,...,N\}$$σ

^{2}is an algorithm-specific variance, and can be tuned to change the extent of parameter-space exploration due to mutation.

**Cloning.**

To ensure that the maximum fitness of the population does not decrease, we copied (cloned) the top few chromosomes with the highest fitness in each generation and inserted them into the next generation.

## 4. Genetic optimization results

#### 4.1. Optimizing planar photonic crystal cavities

### 4.1.1. Q-factor maximization

One problem of interest in PC design is the inverse problem, where one tries to find a dielectric structure to confine a given (target) electromagnetic mode [10]. Here we consider the inverse design problem of optimizing a linear-defect cavity in a planar photonic crystal cavity. The *Q*-factor is a common figure of merit measuring how well a cavity can confine a given mode, and can be approximated (assuming no material absorption) by the following expression:

where *Q*
_{∣∣} represents the *Q*-factor in the direction parallel to the slab, and *Q*
_{⊥} represents the *Q*-factor perpendicular to the slab. *Q*
_{⊥} is usually the limiting factor for *Q _{total}*. As was shown previously [10, 20], the vertical mode confinement, which occurs through total internal reflection (TIR), can be improved if the mode has minimal k-space components inside the light cone.

In the subsequent sections, we report the results for using the GA to minimize the light cone radiation of such cavities. We used one-dimensional photonic crystals as approximations to these cavities [21], and simulated them using the standard Transfer Matrix method for the E-field [22]. The reflectance spectrum of each cavity was obtained using the Transfer Matrix method, and we used a heuristic peak-finding algorithm to automatically search the spectrum for sharp resonance peaks. The resulting resonant modes were then evaluated according to the chosen fitness function (which differed depending on our optimization objective), and the maximum fitness found from all the resonant modes was assigned as the fitness for the particular cavity.

### 4.1.2. Matching to a target function

In [10] it was noted that minimization of light cone radiation could be performed via mode-matching to a target function which already possessed such a property. We therefore used a fitness function that was equal (up to a normalizing factor) to the reciprocal of the mean-squared difference between our simulated mode and a target mode (see Eq (6)). For this simulation, our chromosome encoded the thicknesses of the dielectric slabs in the structure, and was a vector of length 10. We used 100 chromosomes in each generation and allowed them to evolve for 80 generations. Since each dielectric slab was “stacked” above the previous slab, the distances between the centers of the slabs were implicitly coded within our chromosome. We used an alternating piecewise constant dielectric distribution, which varied between *n _{low}* = 1 and

*n*=3.17. The value of

_{high}*n*is an effective refractive index that approximates a true refractive index of a finite-width structure, as will be shown in Section 4.1.4.

_{high}We used target modes that were sinusoidal functions multiplied by *sinc* and *sinc-squared* envelope respectively, which have square and triangular Fourier Transform patterns with no components inside the light cone. Such target modes can have theoretically no radiation in the direction perpendicular to the slab and are therefore ideal candidates as target functions. The results, shown in Fig. 1, clearly feature a suppression of k-vector components at low spatial-frequencies. Matching using the the *sinc-squared* envelope target function produced better results. From the k-space plots, the GA evidently had difficulty matching the sharp edges for the *sinc*-envelope target mode.

### 4.1.3. Direct minimization of light cone radiation

In the preceding subsection, we observed that when we formulated our objective as a matching problem, in the case of the *sinc*-envelope, the GA sacrificed the desired low spatial-frequency suppression in an effort to match the overall shape of the function. The preceding formulation therefore poses an implicit constraint on our optimization. By reformulating the optimization problem, we were able to effectively remove this constraint, and obtain a better result.

Our reformulation directly minimized the k-vector components in the light cone, by minimizing the integrated square-magnitude of the simulated E-field mode in k-space inside the light cone. The fitness function used here is given as in Eq (7), where *V* represents the set of k-vectors within the light cone, and *F*(*k*) is the Fourier Transform of the field.

The final, evolved structure, together with the corresponding real-space and k-space mode profiles are shown in Fig 2. The k-space mode profile features a strong suppression of radiation at low frequencies, to a greater extent as compared to the optimized fields from the preceding simulations. By relaxing our constraint and performing a direct optimization, our GA has designed a structure that achieves better light cone suppression than before. Our direct optimization paradigm has exploited the extreme generality of the GA, which simply requires that a fitness function be defined, with little further constraint thereafter.

### 4.1.4. 2D design verification

In order to verify our design, we used a 2-dimensional Finite-Difference-Time-Domain (FDTD) simulation to compare our GA-optimized design against a standard uniform quarter-wave-stack cavity. We used the standard effective index approximation [23] for uniform slab waveguides to translate our 1-dimensional design (*effective*
*n _{high}* = 3.17) to a 2-dimensional finite-width design (

*true*

*n*= 3.30). Our 2-dimensional cavity therefore has a “slab-waveguide” structure, with a periodic modulation of the dielectric. The uniform cavity comprises a dielectric spacer of thickness λ/2

_{high}*n*with 9 quarter-wave stacks placed symmetrically on each side of the spacer. The number of stacks were chosen so that both the uniform cavity and our GA-optimized cavity have the same number of stacks on each side of the central spacer.

_{high}To isolate the out-of-plane *Q*-factors, we progressively increased the number of quarter-wave-stacks for both the GA-optimized cavity and the uniform cavity, and recorded the overall *Q*-factors. In general, as the number of additional stacks increases, the in-plane confinement of the mode improves, and the overall *Q*-factor converges to the out-of-plane *Q*-factor. There was no significant increase in overall *Q*-factor for either design after adding four additional stacks. The out-of-plane *Q*-factor for the uniform cavity design converged to a value of 85, while the out-of-plane *Q*-factor for our GA-optimized design converged to 6510, surpassing the uniform cavity design by nearly two orders of magnitude. The resulting E-field amplitudes of the resonant modes for both designs are shown in Fig 3. While the uniform cavity displays excellent in-plane modal confinement, it simultaneously exhibits large out-of-plane losses. The authors of Ref [24] noted that this effect was particularly significant for structures with high index ratios, as we have here. Conversely, our GA-optimized cavity effectively trades in-plane modal confinement for a better out-of-plane confinement, leading to a larger overall *Q*-factor.

#### 4.2. Maximal gap at any k-vector point

Moving on to the more general case of 2D photonic crystals, we show the results of simulations for maximizing the TE bandgap at any point in k-space for a 2-Dimensional PC structure with a triangular lattice of air holes. This could be useful for PC design applications where the target mode to be confined is centered around a particular point in k-space [10]. By maximizing the bandgap at that k-space point, we would effectively design a better mirror for a mode resonating along this k-space direction.

We used a supercell which was three periods wide in each dimension (see Fig. 4). We then varied the radii of the nine holes in total, and we encoded the chromosome as a vector of these nine holes. The position of the hole centers were held constant, each spaced apart by one period of the standard triangular lattice. We constrained our search space by restricting each hole’s radius to be less than half a period of the single-hole lattice. This was done to prevent holes from overlapping with each other. A refractive index of *n _{high}*=3.45 was used for the dielectric material, and unity for the refractive index of air. We used a population size of 60 chromosomes for each generation, and allowed the optimization to run for a total of 100 generations. We applied Gaussian mutation to 23 of the 60 chromosomes in each generation, with σ = 0.45

*a*.

### 4.2.1. Maximizing the K-point gap

To evaluate the fitness of each chromosome, we used the eigensolver in Ref [19] to calculate the frequencies of the first 10 bands at the K-point in the reciprocal lattice. Next, we extracted the largest bandgap (calculated as the gap-to-midgap ratio) from these 10 bands. We then scaled the calculated ratio exponentially to tune the selection pressure of the optimization. Figure 5 shows the variation of the gap-to-midgap ratio of our structures as the algorithm progressed.

Our Genetic Algorithm performs as expected, and we get a general increase of fitness as the algorithm progresses. All the four runs do not show any significant increase in fitness after Generation 80, at which point they have maximum fitnesses (i.e. ratio of their bandgap to midgap value) of around 72%. In comparison, a standard uniform triangular lattice, with only one hole per unit cell, has a maximum K-point gap of 53%, achieved when r/a = 0.4445 (maximized by Brent’s algorithm [25]). All the optimized structures have similar dielectric structures and band diagrams; the structures in each of the four runs only differ from each other by a translational shift. The dielectric structures and a sample band diagram are shown in Figure 6. The predicted dielectric structures have, on average, 3 holes per unit cell which have maximum radius (0.5*a*) and 6 holes with negligible radii.

Due to the expansion of the 3 holes and concurrent shrinking of the remaining 6 holes, the superlattice also exhibits discrete translational symmetry in a uniform (albeit *larger*) triangular lattice with periodicity 3*a*. This verifies that the optimal structure for maximizing the K-point gap is still a triangular lattice.

### 4.2.2. Maximizing the M-point gap

To confirm our results, we ran the GA again, this time optimizing the gap at the M-point of the band diagram. Figure 7 shows the results of the optimization. The GA-optimized structures contain, on average, 6 holes per unit cell with large radii close to the upper bound (0.5*a*), and 3 holes with negligible radii. When tiled, the structures appear similar to dielectric waveguides. The M-point gap of the GA-designed structures were about 64%, surpassing the best M-point gap of 55% possessed by a triangular lattice of r/a = 0.4045 (maximized by Brent’s algorithm). Furthermore, we also observe that maximization of the M-point gap comes at the expense of decreasing the K-point gap. These results make intuitive sense, since the waveguide-like structures have a large index contrast in the vertical direction, and therefore can be expected to be vertically well-confined in the M direction, corresponding to a large gap at the M-point.

#### 4.3. Optimal dual PC structures

As a more complex example, let us consider two similar PC designs, (1) a triangular lattice of air holes in a dielectric slab, and (2) a triangular lattice of dielectric rods in air. Structure (1) possesses a bandgap for TE light, but no bandgap for TM light, while structure (2) possesses a bandgap for TM light, but not for TE light.

Our objective is to use the Genetic Algorithm to find a PC design in which the TE eigenmode for structure (1) and the TM eigenmode for structure (2) are most similar. Maxwell’s equations can be cast as eigenproblems for the electric or magnetic fields, and our approach could be potentially useful in future PC design, because solving the inverse problem is analytically simpler (at least intuitively) for the eigenproblem involving the E-field, or for translating TE-like solutions to TM-like solutions.

We used a 3×3 supercell for the optimization, and we minimized the mean-square difference of the z-components of the electric and magnetic fields of the dual structures (oriented along the holes/rods) at the K-point of the band diagram. We encoded the chromosomes as in section 4.2, as a vector of length 9, bounded in the same way as well. However, in this case, each component of the chromosome vector represents *both* the hole radii in structure (1) and the rod radii in structure (2). This ensures that both structures are indeed dual to each other, having an identical geometric arrangement, but with flipped dielectric distributions. We recognize *a priori* that a trivial solution, which we wish to avoid, is a structure that has a uniform refractive index (either dielectric or air) throughout, and so we prevent the genetic algorithm from obtaining this by restricting our mutation to only a Gaussian mutation (see Eq. 4). This preferentially searches the locality of points, and is a necessary trade-off for obtaining a reasonable solution. This illustrates the versatility of the Genetic approach - the extent of the search can be easily modified by a simple change of algorithm parameters. Fig. 8 shows the optimal dual structures with the corresponding simulated fields. The GA clearly converges onto structures that exhibit close similarity between TE and TM confinement. The higher frequency bands begin to deviate because the larger extent of modal variation over the unit cell makes it harder to find a good match.

## 5. Conclusion

We have shown that our Genetic Algorithm is able to effectively optimize PC designs to meet specific design criteria. Specifically, we applied the GA to three particular problems. In 1D cavity simulations, the GA improved vertical cavity confinement by almost two orders of magnitude compared to standard equal-index-spacing designs. We also applied the GA to a 2D triangular lattice to maximize the bandgap at the K-point and M-point of the band diagram. Finally, we use the GA to design symmetric 2D triangular lattice structures that support dual TE and TM modes. Furthermore, by our choice of encoding, we could easily impose constraints upon the design space to ensure that every design searched by the algorithm could be realistically fabricated. Between different optimizations, we only need to change the “fitness function”, which measures how closely a given structure complies with our design criteria. Our Genetic Algorithm is therefore highly robust and can be easily modified to optimize any user-defined objective function.

## Acknowledgements

This work has been supported by the MURI Center for photonic quantum information systems (ARO/DTO Program DAAD19-03-0199), and NSF Grant No. CCF-0507295.

## References and links

**1. **S. John, “Strong localization of photons in certain disordered dielectric superlattices,” Phys. Rev. Lett. **58**, 2486–2489 (1987). [CrossRef] [PubMed]

**2. **E. Yablonovitch, “Inhibited spontaneous emission in solid-state physics and electronics,” Phys. Rev. Lett. **58**, 2059–2062 (1987). [CrossRef] [PubMed]

**3. **D. Englund, D. Fattal, E. Waks, G. Solomon, B. Zhang, T. Nakaoka, Y. Arakawa, Y. Yamamoto, and J. Vučković, “Controlling the spontaneous emission rate of single quantum dots in a two-dimensional photonic crystal,” Phys. Rev. Lett. **95** (2005). [CrossRef] [PubMed]

**4. **M. Boroditsky, R. Vrijen, T. Krauss, R. Coccioli, R. Bhat, and E. Yablonovitch, “Control of spontaneous emission in photonic crystals,” Proceedings of SPIE - The International Society for Optical Engineering **3621**, 190–197 (1999).

**5. **H. Altug and J. Vučković, “Experimental demonstration of the slow group velocity of light in two-dimensional coupled photonic crystal microcavity arrays,“ Appl. Phys. Lett. **86** (2005). [CrossRef]

**6. **Y. A. Vlasov, M. O’Boyle, H. F. Hamann, and S. J. McNab, “Active control of slow light on a chip with photonic crystal waveguides,” Nature **438**, 65–69 (2005). [CrossRef] [PubMed]

**7. **H. Altug and J. Vučković, “Photonic crystal nanocavity array laser,“ Opt. Express **13**, 8819–8828 (2005). [CrossRef] [PubMed]

**8. **B.-S. Song, S. Noda, T. Asano, and Y. Akahane, ”Ultra-high-Q photonic double-heterostructure nanocavity,“ck Nat. Mat. **4**, 207–210 (2005). [CrossRef]

**9. **J. Vučković, M. Lončar, H. Mabuchi, and A. Scherer, “Design of photonic crystal microcavities for cavity QED,” Phys. Rev. E **65** (2001).

**10. **D. Englund, I. Fushman, and J. Vucčković, “General recipe for designing photonic crystal cavities,” Opt. Express **13**, 5961–5975 (2005). [CrossRef] [PubMed]

**11. **D.A.B. Miller, Y. Jiao, and S. Fan, “Demonstration of systematic photonic crystal device design and optimization by low-rank adjustments: an extremely compact mode separator,” Opt. Lett. **30**, 141–143 (2005). [CrossRef] [PubMed]

**12. **S. Preble, H. Lipson, and M. Lipson, “Two-dimensional photonic crystals designed by evolutionary algorithms,” Appl. Phys. Lett. **86** (2005). [CrossRef]

**13. **R. P. Drupp, J. A. Bossard, D. H. Werner, and T. S. Mayer, “Single-layer multiband infrared metallodielectric photonic crystals designed by genetic algorithm optimization,” Appl. Phys. Lett. **86** (2005). [CrossRef]

**14. **J. H. Holland, *Adaptation in Natural and Artificial Systems: An Introductory Analysis with Applications to Biology, Control and Artificial Intelligence*, Univ. of Michigan Press (1975). [PubMed]

**15. **D. E. Goldberg, *Genetic Algorithms in Search, Optimization and Machine Learning*, Addison Wesley (1989).

**16. **L. Davis, *Genetic Algorithms and Simulated Annealing*, Morgan Kaufmann (1987).

**17. **L. Shen, Z. Ye, and S. He, “Design of two-dimensional photonic crystals with large absolute band gaps using a genetic algorithm,” Phys. Rev. B **68** (2003). [CrossRef]

**18. **E. Kerrinckx, L. Bigot, M. Douay, and Y. Quiquempois, “Photonic crystal fiber design by means of a genetic algorithm,“ Opt. Express **12**, 1990–1995 (2004). [CrossRef] [PubMed]

**19. **S. G. Johnson and J. D. Joannopoulos, “Block-iterative frequency-domain methods for maxwell’s equations in a planewave basis,” Opt. Express **8**, 173–190 (2001). [CrossRef] [PubMed]

**20. **Y. Akahane, T. Asano, B.-S. Song, and S. Noda, “High-*Q* photonic nanocavity in a two-dimensional photonic crystal,” Nature **425**, 944–947 (2003). [CrossRef] [PubMed]

**21. **P. Lalanne, S. Mias, and J. Hugonin, “Two physical mechanisms for boosting the quality factor to cavity volume ratio of photonic crystal microcavities,” Opt. Express **12**, 458–467 (2004). [CrossRef] [PubMed]

**22. **A. Yariv and P. Yeh, *Optical Waves in Crystals: Propagation and Control of Laser Radiation*, John Wiley and Sons Inc (2002).

**23. **L. Coldren and S. Corzine, *Diode Lasers and Photonic Integrated Circuits*, John Wiley and Sons Inc (1995).

**24. **J. Vučković, M. Pelton, A. Scherer, and Y. Yamamoto, “Optimization of three-dimensional micropost microcav-ities for cavity quantum electrodynamics,” Phys. Rev. A **66** (2002).

**25. **R. Brent, *Algorithms for Minimization Without Derivatives*. Prentice-Hall (1973).