OSA's Digital Library

Applied Optics

Applied Optics

APPLICATIONS-CENTERED RESEARCH IN OPTICS

  • Editor: Joseph N. Mait
  • Vol. 48, Iss. 21 — Jul. 20, 2009
  • pp: 4177–4190
« Show journal navigation

Detailed electromagnetic simulation for the structural color of butterfly wings

R. Todd Lee and Glenn S. Smith  »View Author Affiliations


Applied Optics, Vol. 48, Issue 21, pp. 4177-4190 (2009)
http://dx.doi.org/10.1364/AO.48.004177


View Full Text Article

Acrobat PDF (2883 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

Many species of butterflies exhibit interesting optical phenomena due to structural color. The physical reason for this color is subwavelength features on the surface of a single scale. The exposed surface of a scale is covered with a ridge structure. The fully three-dimensional, periodic, finite-difference time- domain method is used to create a detailed electromagnetic model of a generic ridge. A novel method for presenting the three-dimensional observed color pattern is developed. Using these tools, the change in color that is a result of varying individual features of the scale is explored. Computational models are developed that are similar to three butterflies: Morpho rhetenor, Troides magellanus, and Ancyluris meliboeus.

© 2009 Optical Society of America

1. Introduction

Many species of butterfly display vibrant iridescent coloration that is widely believed to be caused by optical interference. This type of coloration is called structural color (see [1

1. S. Berthier, Iridescences: The Physical Colors of Insects (Springer, 2007).

] for an overview of the field). Structural color can have several interesting features: The hue often changes with viewing angle, the color is often very intense and highly saturated, and the color can be changed or even eliminated by immersing the insect in a liquid with an appropriate index of refraction.

Scientific interest in structural color of insects has persisted for well over a century. Newton, Rayleigh, and Michelson all contributed theories to the field. The major breakthroughs in understanding occurred with the invention of the electron microscope [2

2. T. F. Anderson and A. G. Richards, “An electron microscope study of some structural colors of insects,” J. Appl. Phys. 13, 748–758 (1942). [CrossRef]

]. For the first time, it became possible to observe the subwavelength features that are the source of structural color. Talented biologists, notably Ghiradella, published extensive studies of the structural features [3

3. H. Ghiradella, “Structure of iridescent lepidopteran scales: variations on several themes,” Ann. Entomol. Soc. Am. 77, 637–645 (1984).

, 4

4. H. Ghiradella, “Hairs, bristles, and scales,” in Microscopic Anatomy of Invertebrates, W. H. Frederick and L. Michael, eds. (Wiley-Liss, 1998), Vol. 11A, pp. 257–287.

]. Her studies, among others, were used to create a classification scheme for insects based on the underlying optical principles of the color [5

5. P. Vukusic, J. R. Sambles, and H. Ghiradella, “Optical classification of microstructure in butterfly wing-scales,” Photonics Sci. News 6, 61–68 (2000).

].

A typical butterfly exhibiting structural color is the male Morpho rhetenor, shown in Fig. 1. This species is a very intense blue that fades to a deep violet when viewed near grazing. The wings are covered with a large number of small scales that are easily visible using a low-power optical microscope. As shown in Fig. 1b, it is actually the scales that produce the color. The scales, in turn, are covered with a series of ridges that are less than 1μm apart, which are shown in Fig. 1c. Each ridge is a treelike structure with 4–12 branches on both sides. The structure is made of cuticle with an approximate complex index of refraction n˜=1.56+i0.06 [6

6. P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wooten, “Quantified interference and diffraction in single Morpho butterfly scales,” Proc. R. Soc. London Ser. B 266, 1403–1411 (1999). [CrossRef]

].

In some species of butterfly, there are two types of scales. The first type, known as a ground scale, is similar to the scales on the M. rhetenor. In addition, a second type of scale, known as a cover scale, is layered over the ground scales. These scales generally have fewer ridges than the ground scales and are semitransparent. They are thought to diffuse the scattering from the ground scales [3

3. H. Ghiradella, “Structure of iridescent lepidopteran scales: variations on several themes,” Ann. Entomol. Soc. Am. 77, 637–645 (1984).

]. This work primarily concentrates on ground scales, although a simplified cover scale is considered.

In this paper we use one of the most powerful numerical techniques for solving electromagnetic and optical problems, the finite-difference time-domain (FDTD) method [7

7. A. Taflove and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech, 2005).

], to study the scattering of light by the ridge structure described above. We perform a fully three-dimensional analysis that makes use of the most recent refinements to the method, which include a new approach for handling periodic structures [8

8. R. T. Lee and G. S. Smith, “An alternative approach for implementing periodic boundary conditions in the FDTD method using multiple unit cells,” IEEE Trans. Antennas Propag. 54, 698–705 (2006). [CrossRef]

, 9

9. R. T. Lee, “A novel technique for incorporating periodic boundaries into the FDTD method and the application to the study of structural color of insects,” Ph.D.thesis(Georgia Institute of Technology, 2009).

]. In the past, the FDTD method has been used to study butterfly scales; however, restrictions not introduced here, such as a reduced number of dimensions (two), were employed [10

10. L. Plattner, “Optical properties of the scales of Morpho rhetenor butterflies: theoretical and experimental investigation of the backscattering of light in the visible spectrum,” J. R. Soc. Interface 1, 49–59 (2004). [CrossRef]

, 11

11. S. Banerjee, J. B. Cole, and T. Yatagai, “Colour characterization of a Morpho butterfly wing-scale using a high accuracy nonstandard finite-difference time-domain method,” Micron 38, 97–103 (2007). [CrossRef]

]. In this paper we also present a novel method for displaying the color associated with the three-dimensional scattering pattern for the scale. We vary several of the features of the scale and observe the changes in color that result by using this new method. Our goal is not to create an exact computational model for the M. rhetenor or any other butterfly; rather, it is to study the effect of varying geometrical parameters on a generic scale.

2. Computational Model

2A. Discretized Geometry

By examining electron micrographs such as Fig. 1c, a computational model for a generic butterfly scale can be developed. An example of such a model is shown in Fig. 2a. In this model, the base, which is the relatively featureless cuticle sheet at the bottom of Fig. 1c, is modeled as a dielectric slab. The upper structure, known as a ridge, is a set of dielectric slabs with different dimensions and orientations. First, there is a central stem that supports the structure. On either side of this stem are branchlike arms called lamellae. The lamellae are often tapered, with those closest to the base being the widest. Usually, the lamellae are slightly tilted in the longitudinal direction (y) with respect to the base; however, some species have a very high tilt, as is discussed below. In addition to the lamellae, there is sometimes a second set of layers nearly orthogonal to the first. These layers are smaller than the lamellae and are called microribs. In Fig. 2, the components of the ridge are shown in different colors for ease of identification. In the model, they all have the same electrical properties: electrical permittivity ϵ=2.43ϵo, magnetic permeability μ=μo, and electrical conductivity σ=6.4×103S/m. The permittivity and conductivity were obtained from the aforementioned complex index of refraction, n˜, at the wavelength λ=488nm.

The structure in Fig. 2a is doubly periodic. The periodic cell, indicated by the black frame, is repeated three times in each lateral direction in this figure. The periodic cell is examined in more detail in Fig. 3. In this figure, the width of the periodic cell, xp, is the period of the ridges, and the length, yp, is determined by the spacing between the lamellae, g, the thickness of the lamellae, h, and the tilt angle of the lamellae with respect to the base, α, according to
yp=(g+h)/sinα.
(1)
If present, the thickness and spacing of the microribs, f and e, respectively, are constrained so that an integer number of microribs lie within a single periodic cell.

Figures 3a, 3b, 3c show three variations on the generic model. The number of lamellae in all cases is N. In (a) the lamellae are all of the same width, a, and those on the left and right sides of the stem are aligned. In (b) some of the lamellae are tapered, so the maximum width is a and the minimum width is d. Finally, in (c) the lamellae on the left and right sides of the stem are offset vertically instead of aligned. Notice that there is also a dielectric layer above the ridges in (c). It is intended to be a highly simplified model for the aforementioned cover scale.

In the actual scales, there is also structure between the ridges and the base that supports the ridges. Our calculations show that this structure has little effect on the color of the scattered light, so it is omitted in this model.

2B. Finite-Difference Time-Domain Method

At this time, the FDTD method has been developed to a point at which many of the features are quite complicated. Our purpose here is not to cover the details of these features; rather, it is to give a simple explanation of the basic algorithm and an overview of the special features used for the current problem.

The FDTD method is used to directly solve Maxwell’s equations in the time domain [7

7. A. Taflove and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech, 2005).

]. In time domain, vector form, the equations are
×=μHt,
(2)
×H=σ+ϵt,
(3)
in which E is the electric field and H is the magnetic field. Throughout the paper, a quantity in the time domain is indicated by script font, such as E, and a quantity in the frequency domain is indicated by a Roman font, such as E. The volume to be simulated is divided into a large number of electrically small blocks (Δx×Δy×Δz), as shown in Fig. 2b. Typically, the dimensions for the blocks are less than λmin/20, where λmin is the shortest wavelength in the simulation (this is determined by the excitation). The electromagnetic field is evaluated at discrete points in space (specified by indices i, j, k) and time (index n). For each block, the vector components of the electromagnetic field are stored. The components are offset in space, as shown in Fig. 2c. They are also offset in time so that the electric field components are known at time t=nΔt, and the magnetic field components are known at t=(n+1/2)Δt. This particular arrangement of field components is known as a Yee cell [12

12. K. S. Yee, “Numerical solution of initial boundary value prob lems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966). [CrossRef]

]. The standard FDTD notation for a field component is
Hz|i,j,kn+1/2=Hz(iΔx,jΔy,kΔz,(n+1/2)Δt).
(4)

The FDTD algorithm uses a procedure called “marching in time” to calculate future values of the electromagnetic field from current and past values. For example, the value of Hz at the future time t=(n+1/2)Δt is determined from the value of Hz at the past time t=(n1/2)Δt and the values of Ex and Ey at the current time t=nΔt using the update equation
Hz|i,j,kn+1/2=Hz|i,j,kn1/2+ΔtμΔy(Ex|i,j+1/2,knEx|i,j1/2,kn)ΔtμΔx(Ey|i+1/2,j,knEy|i1/2,j,kn).
(5)
Similar equations are used for the remaining field components. The algorithm is documented in detail in [7

7. A. Taflove and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech, 2005).

].

Figure 4 is a schematic drawing for the FDTD computational domain used for our problem, one periodic cell. A field must be introduced into this cell before the algorithm in Eq. (5) can be used to advanced the solution in time. For this problem, the field is an incident plane wave, and it is produced by the “plane-wave injector” near the top of the cell. Here we use the recently developed “nearly perfect” plane-wave injector described in [13

13. J. B. Schneider, “Plane waves in FDTD simulations and a nearly perfect total-field/scattered-field boundary,” IEEE Trans. Antennas Propag. 52, 3280–3287 (2004). [CrossRef]

]. The injector surface is designed to produce the incident field below, toward the structure, and no field above. As a result, the injector surface naturally divides the domain into two regions: the region below in which the incident field and the scattered field (the field produced by the interaction of the incident field with the structure) exist, and the region above in which only the scattered field exists. These two regions are referred to as the “total field region” and “scattered field region,” respectively.

The incident field produced by the injector is a linearly polarized plane wave with an electric field of the form
Ei(r,t)=G(tk^·r/c)e^,
(6)
in which k^ is a unit vector in the direction of propagation, and e^ is the polarization vector. The time-varying scalar function G must have sufficient bandwidth to cover the optical frequencies (wavelengths) of interest. For this work, it is a sinusoid modulated by a Gaussian pulse.

2C. Modeling of Natural Light

We will consider butterfly scales exposed to natural light (daylight). Natural light is incoherent light that is unpolarized. Here it is modeled by making two separate calculations (1 and 2) in which the states of linear polarization of the incident plane waves, represented by Eq. (6), are orthogonal: e^1·e^2=0. All subsequent calculations based on these fields, such as reflection and transmission coefficients, involve the power, and they are obtained by averaging the appropriate powers for the two states.

3. Effect of Randomness in the Geometry

A periodic structure, such as the one shown Fig. 2, is a diffraction grating. When it is infinite in lateral extent, the scattered radiation is a discrete set of plane waves, the directions of which are determined by a grating equation. For a large grating (many periods), the scattered radiation consists of a number of very narrow lobes (grating lobes), with their directions again given by the grating equation. Observations of the scattered light from butterfly scales and wings do not show these narrow lobes. Using their measurements, Kinoshita and Yoshioka have determined that the absence of these lobes is probably due to irregularities in the structure, namely, random positioning of the ridges on a scale [15

15. S. Kinoshita and S. Yoshioka, “Structural colors in nature: The role of regularity and irregularity in the structure,” Chem. Phys. Chem. 6, 1442–1459 (2005). [CrossRef] [PubMed]

]. As a result of this randomness, the scattered fields from the ridges add incoherently, rather than coherently as in a perfect diffraction grating. The assumption of randomness, which we adopt, implies that the scattering from a whole scale can be determined from that of a single ridge. Below we describe numerical experiments that were performed to better understand this assumption and how to implement it in the FDTD method.

To reduce the computational burden (storage and run time), the numerical experiments were performed using a simplified version of the model shown in Fig. 2. The model was two dimensional; that is, it was invariant in the y direction (α=0, yp) with a cross section similar to the one shown in Fig. 3b with the top four lamellae tapered, but without the microribs. The other dimensions were N=8, c=60nm, a=620nm, d=140nm, h=60nm, g=140nm, and xp=750nm.

First, a random array was simulated using n=50 ridges. The vertical positions of the ridges, the displacements Δh in Fig. 5a, were varied randomly with a triangular probability density, such that 275nmΔh275nm. The far-zone or Fraunhofer scattered field (limkr) of the array was determined using a form of Huygens’ principle, which states that the far-zone scattered field can be determined from an integral involving the near-zone scattered field on a closed surface surrounding the array, such as the box marked “transform surface” in Fig. 5a. In the literature for the FDTD method, this procedure is called a near-to-far-field transformation. Now for a large array, the right and left sides of the box make relatively small contributions to the integral and can be ignored in the calculation. In addition, the total field below the array is approximately zero; all of the energy in the incident field is scattered or absorbed before reaching this surface. Thus the scattered field on the bottom surface of the box is approximately the negative of the incident field, which is unaffected by the randomness in the array. This coherent field on the bottom of the box produces a beam of radiation below the box (approximately the negative of the incident field) and little radiation above the box. From these observations, we conclude that it is mainly the near-zone scattered field on the top surface of the box that determines the far-zone scattered field above the box (the backscattered field). In the calculations that follow, only this surface will be used.

To obtain a statistical representation, the far-zone scattered field from the array of 50 randomly placed ridges was calculated 100 times using different random displacements for each case. The ensemble average of the 100 cases was then computed. This average can be thought of as what one sees when viewing a section of a wing containing many scales. The graphs in Fig. 6a are the results of those calculations, which are far-zone scattering patterns that show the time-average scattered power per unit area as a function of the direction θ in the x–z plane. Notice that the scattering is largest at the wavelengths that represent blue–violet light. The small “noiselike” variations along each curve are due to the random aspect of the calculation.

Next, we show how we can obtain results similar to those for the random array, shown in Fig. 5a, using the periodic FDTD method. For determining the far-zone scattered field above the array, the surface shown in Fig. 5b is equivalent to the top surface of the box in Fig. 5a. We ignore the small vertical portions of this surface; that is, we only consider the horizontal surfaces marked S1,S2,...,Sn in our calculation. Because of the random displacement of the ridges, the scattering from these surfaces is incoherent, and in a statistical sense, it is equivalent to n times the scattering from one of the surfaces. Thus, knowledge of the near-zone scattered field on the surface immediately above one ridge can be used to obtain the far-zone scattered field from the random array of ridges. We approximate this field by the field from the periodic FDTD calculation based on the unit cell shown in Fig. 4, the field on the “transform surface” in the scattered field region. We note that because of the periodic boundary conditions in this calculation, the field within the unit cell includes the influence of all of the other ridges in the array.

The results based on the periodic FDTD calculation are shown in Fig. 6b, and they are to be compared to those for the random array in Fig. 6a. The scattering patterns for the two methods are very similar (the former is a smoothed version of the latter), with the shapes of the curves and relative values for the different wavelengths being essentially the same.

The power reflection coefficient is a measure of the total power scattered per ridge of the scale. It is obtained by integrating the normal component of the Poynting vector over the transform surface, then dividing by the total power incident on this surface. We note that, for the random array, the transform surface is the area above the 50 ridges that make up one element in the ensemble average, whereas for the periodic element, the transform surface is the area above the ridge (unit cell). The power reflection coefficients for the two methods are compared in Fig. 6c. As for the scattering patterns, we see that the results from the two methods are very similar.

The results presented in this section show that FDTD calculations based on the single periodic cell can be used to estimate the scattering, in an average sense, from a scale containing a large number of irregularly displaced ridges.

4. Interpretation of Scattering as Observed Color

In Section 3, we described the use of the periodic FDTD method to calculate the scattered field from a scale with irregularly positioned ridges. In this section, we explain the procedure used to obtain the observed color of a butterfly model from the results of the FDTD simulations. The goal is to calculate the observed color as a function of observer location when the structure is illuminated by an unpolarized light source located infinitely far away.

Referring to the flowchart in Fig. 7, the first step is to specify the incident field. Two simulations are performed using plane-wave excitations with orthogonal states of polarization. The scattered field is recorded on the transform surface above the structure. Although the computational method is a time- domain method, frequency-domain quantities are desired. Therefore the scattered field is converted to the frequency domain using an on-the-fly discrete Fourier transform (DFT) [7

7. A. Taflove and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech, 2005).

]. The frequencies used correspond to 81 wavelengths spaced in 5nm increments between λ=380nm and λ=780nm. This set of wavelengths covers the visible spectrum.

After the simulations are completed, the far-zone scattered field is calculated using the previously described frequency-domain, near-to-far-field transform (Huygens’ principle). From the far-zone field, the radiant intensity, I, of the scattered light is found. In addition, the radiant flux of the incident light, Ii, is calculated. Note that both quantities are functions of the wavelength, but the scattered radiant intensity is also a function of the direction. The reflection factor, R(λ,θ,ϕ), is calculated as the ratio of the scattered radiant intensity to the incident radiant flux [16

16. “Standard terminology of appearance,” ASTM E 284 08 (ASTM International, 2008).

].

At each observation direction, the reflection factor is converted into tristimulus values (X,Y,Z), for example,
X(θ,ϕ)=KcλD65(λ)R(λ,θ,ϕ)x¯(λ)Δλ,
(7)
in which D65 is the spectrum of a standard reference white light source with an approximate color temperature of 6500K, and x¯ is a color matching function. For this work, the 1931 CIE 2° standard observer color matching functions are used [17

17. R. S. Berns, Billmeyer and Saltzman’s Principles of Color Technology, 3rd ed. (Wiley, 2000).

]. By using these functions, the tristimulus value Y corresponds to the luminance factor. The constant Kc is set so that a perfect diffuser, that is, R(λ,θ,ϕ)=constant, is white. Although the (X,Y,Z) values uniquely specify a color, it is more convenient to convert to standard red-green-blue (sRGB) values that can be displayed on a standard computer monitor [16

16. “Standard terminology of appearance,” ASTM E 284 08 (ASTM International, 2008).

, 18

18. “Standard practice for computing the colors of objects by using the CIE system,” ASTM E 308 06 (ASTM International, 2006).

].

The observed color is plotted as a three-dimensional pattern; the procedure is shown schematically in Fig. 8. The radial vector from the origin O is proportional to the luminance factor Y. The small patch of surface (angular size Δθ, Δϕ) at the end of this vector is the color that would be seen by an observer viewing the structure from that direction (θ,ϕ). By examining such patterns, one can see how the luminance and color of the scattered light changes with the direction of observation.

5. Results

5A. Models for Morpho Butterflies

In this subsection, we examine the scattering from structures that are similar to a Morpho rhetenor scale (once again, the computational model is inspired by the specific insect, but it is not meant to be an exact model). The dimensions of the structure are based on the observations described in [19

19. R. A. Potyrailo, H. Ghiradella, A. Vertiatchikh, K. Dovidenko, J. R. Cournoyer, and E. Olson, “Morpho butterfly wing scales demonstrate highly selective vapour response,” Nat. Photon. 1, 123–128 (2007). [CrossRef]

] and are summarized in the first column of Table 1.

Figure 9 shows the observed color for normally incident light (θi=0°) as a function of the location of the observer (θ,ϕ). Each subplot corresponds to a change in the structure. Notice the small insets in Figs. 9a, 9d; they are simplified models that show the orientation of the ridge. In this figure, and in all similar figures that follow, logarithmic scaling is used; that is, the length of the radial vector is proportional to log10(Y/Ymax), and the range displayed is two decades. Here, Ymax is the maximum value for each model for a given genus, and it is given for each subplot in the caption of the figure. Using these values, we can compare the relative luminance values for the subplots.

Figure 9a is for the simplest model; it is the one in Fig. 3a with a base and without microribs. The eight lamellae on the left and right sides of the stem are aligned, and all have the same width. Notice that the color is a highly saturated blue and is concentrated in a small angular sector, a beam. This beam is tilted at approximately the angle θ=2α=20°. That is, the tilt of the lamellae with respect to the plane of the wing (base) causes a rotation in the beam that is approximately twice the tilt angle. This rotation is the same as one would expect for specular reflection from a tilted dielectric slab.

To provide a quantitative measure for the shape of the beam in Fig. 9a, the scattered power is examined for two planes. Figure 10a is for the plane of the ridge (y–z plane), while Fig. 10b is for a plane perpendicular to the y–z plane and passing through the peak of the beam. Curves are given for four wave lengths in the visible range. The aforementioned tilt of the beam to θ20° is evident in Fig. 10a. For all wavelengths, the width of the beam is seen to be significantly greater in the perpendicular plane, Fig. 10b. This behavior is consistent with measured results for Morpho butterflies [6

6. P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wooten, “Quantified interference and diffraction in single Morpho butterfly scales,” Proc. R. Soc. London Ser. B 266, 1403–1411 (1999). [CrossRef]

, 20

20. S. Berthier, E. Charron, and J. Boulenguez, “Morphological structure and optical properties of the wings of Morphidae,” Insect Sci. 13, 145–158 (2006). [CrossRef]

].

The power reflection coefficient versus wavelength for this case is shown in Fig. 11a. There is a definite peak at approximately λ470nm, which is in the blue range of the visible spectrum. The reflection coefficient is similar to the measurements of Berthier ([1

1. S. Berthier, Iridescences: The Physical Colors of Insects (Springer, 2007).

], Fig. 8.20).

For Fig. 9b, the structure is identical to that for Fig. 9a, except the top four lamellae are tapered in width, as in Fig. 3b. The pattern of the scattered light is similar for both cases. However, as the power reflection coefficient in Fig. 11a shows, the intensity of the light is lower with the taper. This is also evident from the luminance scale factors for Figs. 9a, 9b: Ymax=1.0 versus Ymax=0.30.

For Fig. 9c, the model is changed once again by adding a cover scale (dielectric slab) above the ridge. Recall that the cover scale is prevalent in Morpho-type butterflies but not the M. rhetenor [6

6. P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wooten, “Quantified interference and diffraction in single Morpho butterfly scales,” Proc. R. Soc. London Ser. B 266, 1403–1411 (1999). [CrossRef]

]. Two beams can now be identified in the pattern, a beam due to specular reflection from the ridge (θ=2α=20°), which is similar in shape and hue to the one in Fig. 9b, and a beam due to specular reflection from the diffusing layer (θ=0°), which is of different hue and less saturated. The combination produces a pattern that is broader than in Fig. 9b, which is consistent with the concept of a diffusing layer. The decrease in saturation is consistent with observations of Morpho butterflies with cover scales, such as M. didius, for which the wings are less glossy than the M. rhetenor, and the observed color is whitish-blue [21

21. S. Yoshioka and S. Kinoshita, “Wavelength-selective and anisotropic light-diffusing scale on the wing of the Morpho butterfly,” Proc. R. Soc. London Ser. B 271, 581–587 (2004). [CrossRef]

].

For Fig. 9e, the top four lamellae are tapered in width. The taper is seen to shift the hue slightly toward violet. The effect of the taper on the intensity is more dramatic and similar to that seen for aligned lamellae. It causes a significant decrease, as can be seen in the power reflection coefficient in Fig. 11b and in the scale factors for the luminance: Ymax=0.05 versus Ymax=0.21.

Finally, the effect of a cover scale (dielectric slab) is shown in Fig. 9f. The specular reflection (beam) from the cover scale is significantly larger than the scattering from the ridge below and produces a fairly unsaturated blue color. Comparing this result with the one for aligned lamellae in Fig. 9c, we can directly see the effect of offsetting the lamellae. The scattering from the cover scale and ridge are comparable when the lamellae are aligned; however, offsetting the lamellae significantly reduces the scattering from the ridge, making the scattering from the cover scale the dominant feature.

Notice that in all cases shown in Figs. 9, there is a change in hue with viewing position (the angles θ and ϕ). This is a clear indication of the iridescence observed for these butterflies.

Next, the effect of adding microribs to the Morpho-like models is examined. We recall that the microribs are the small layers that are perpendicular to the lamellae, as shown in Figs. 2, 3. The microribs are added to the models with aligned and offset lamellae with taper; results for these cases are presented in Figs. 9b, 9e. The results in Fig. 12 show that the addition of the microribs causes very little change in the shape of the patterns and only a small change in the hue for the Morpho-like models. The reason is straightforward: The lamellae are tilted at a slight angle (10°) relative to the base; thus the microribs are nearly perpendicular to the base. Therefore specular scattering from the microribs is mainly reflected into the adjacent structure rather than back at the observer.

5B. Troides magellanus Butterfly

Next, the scattering from a structure similar to the male Troides magellanus butterfly scale is studied. This insect, shown in Fig. 13a, is remarkable due to the interplay of structural and pigmentary color [23

23. C. Lawrence, P. Vukusic, and R. Sambles, “Grazing-incidence iridescence from a butterfly wing,” Appl. Opt. 41, 437–441 (2002). [CrossRef] [PubMed]

, 24

24. J. P. Vigneron, K. Kertesz, Z. Vertesy, M. Rassart, V. Lousse, Z. Balint, and L. P. Biro, “Correlated diffraction and fluorescence in the backscattering iridescence of the male butterfly Troides magellanus (Papilionidae),” Phys. Rev. E 78, 021903 (2008). [CrossRef]

]. When viewed at any angle other than near grazing, the lower wings of the butterfly are yellow. However, when the observer and the light source are collocated in a direction near grazing and pointing toward the outer margin of the wing, the observed color changes with angle to blue, as shown in Fig. 13b.

The parameters used in the model for this butterfly are given in the second column of Table 1; they are based on values given in [23

23. C. Lawrence, P. Vukusic, and R. Sambles, “Grazing-incidence iridescence from a butterfly wing,” Appl. Opt. 41, 437–441 (2002). [CrossRef] [PubMed]

, 24

24. J. P. Vigneron, K. Kertesz, Z. Vertesy, M. Rassart, V. Lousse, Z. Balint, and L. P. Biro, “Correlated diffraction and fluorescence in the backscattering iridescence of the male butterfly Troides magellanus (Papilionidae),” Phys. Rev. E 78, 021903 (2008). [CrossRef]

]. The key feature is that the tilt angle of the lamellae is large, α54°. Recall that the tilt angle was only 10° for the Morpho-like structures. Due to the larger tilt, the period along the y direction, yp, is shorter than for the Morpho models by about a factor of four. Therefore, to keep the length for coherence the same for both butterflies, four periods along the y direction are used for the calculations for this insect. The other parameters for the lamellae are similar to those for the Morpho-like models; therefore it is expected that the observed color will be similar under certain conditions. Specifically, when the direction of the incident light is near normal to the lamellae (approaching grazing to the wing), then the observed color should be blue.

This hypothesis is confirmed by examining the plots in Figs. 14, 15. When the angle of incidence is well away from grazing, the scattering from the structure is very weak. This is the case in Fig. 14a, which is for θi=40°. Notice that the scale factor for the luminance is Ymax=0.05. The color is brown, not yellow as for the actual butterfly, because the yellow color is due to pigmentation, which is not included in the model. As the angle of incidence is increased slightly toward grazing, a lobe that is distinctly blue in color appears, as in Fig. 14b, which is for θi=50°. This blue lobe is even more pronounced when θi=70°, as in Fig. 14c. Notice that the scale factor for the luminance is now Ymax=1.0. As for the Morpho, an intense blue appears when the light is nearly normally incident on the broad surfaces of the lamellae.

There are factors in addition to the direct scattering from the broad surfaces of the lamellae that contribute to the color shown in Fig. 14. As pointed out by Lawrence et al., the exposed ends of the highly tilted lamellae form a grating on the horizontal (x–y) plane [23

23. C. Lawrence, P. Vukusic, and R. Sambles, “Grazing-incidence iridescence from a butterfly wing,” Appl. Opt. 41, 437–441 (2002). [CrossRef] [PubMed]

]. For our dimensions, the period of this grating is pe=(g+h)/sinα=260nm, and there are two lobes in the scattering pattern. One is a specular lobe (zeroth order) that is fairly wavelength independent, and is probably responsible for the dark lobe on the right of all of the plots in Fig. 14. The other lobe (first order) occurs close to the backscattering direction (θ=θi) for wavelengths shorter than λ=2pe=520nm, and it contributes to the blue-violet color.

In Fig. 15, the power reflection coefficient is plotted versus wavelength for two angles of incidence, θi=40° and θi=70°. These results clearly show the increase in scattering that occurs as the angle of incidence approaches grazing. They also clearly indicate the blue color for the structural scattering.

5C. Ancyluris meliboeus Butterfly

The observed color for the model with microribs added is shown in Fig. 16b. The shape of the beam, which is largest near the specular direction (θ=60°), is not significantly changed, but the color is shifted toward longer wavelengths. In addition, a second lobe appears. This lobe is believed to be caused by the scattering from the microribs acting as a second set of lamellae. The shift in color with the addition of the microribs is also indicated by the power reflection coefficients in Fig. 17. The observed color is similar to that in photographs of a specimen taken by Vukusic et al. [25

25. P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wootton, “Limited-view iridescence in the butterfly Ancyluris meliboeus,” Proc. R. Soc. London Ser. B , 269, 7–14 (2002). [CrossRef]

].

6. Conclusion

In this work, a set of detailed electromagnetic simulations for the structural color of butterfly wings was presented. The simulations were performed using the finite-difference time-domain method with periodic boundary conditions. The effect of the irregularity between the ridges of a scale was examined, and a method for using a periodic model to predict the response from an irregular structure was developed. The results from the simulations were interpreted using colorimetric theory in terms of human color vision. A novel method for presenting the three- dimensional color pattern was introduced. The techniques were used to study computational models that are similar to the structures of three butterflies: Morpho rhetenor, Troides magellanus, and Ancyluris meliboeus. The effect of variations in the structure of these models was explored.

This work was supported in part by the John Pippin Chair in Electromagnetics within the School of Electrical and Computer Engineering at the Georgia Institute of Technology.

Table 1. Dimensions for Butterfly Models

table-icon
View This Table
Fig. 1 Images of a male Morpho rhetenor butterfly. (a) Photograph of dorsal side (top) of insect. (b) Detail of wing showing scales as viewed under an optical microscope. (c) Cross section of scale showing individual ridges as viewed using a scanning electron microscope. SEM image courtesy S. Kinoshita [[26], p. 13, Fig. 12].
Fig. 2 Computational model. (a) Portion of a scale showing three ridges. The volume within the black frame is the volume used for the periodic FDTD calculations. (b) Detail of a single ridge showing the individual FDTD cells. (c) Arrangement of the electromagnetic field components within a single FDTD cell.
Fig. 3 Parameters used to describe the geometry of a periodic cell and the coordinates for a point P in space.
Fig. 4 Schematic drawing showing elements within an FDTD periodic cell.
Fig. 5 (a) Geometry for two-dimensional array of ridges with random heights. (b) Equivalent transform surface for obtaining the far-zone scattered field above the array.
Fig. 6 Comparison of results from an ensemble average of randomized finite models and a periodic model. (a) Far-zone scattering pattern for the ensemble average of 100 simulations with 50 randomly spaced ridges per simulation. (b) Far-zone scattering pattern for the periodic model. (c) Comparison of the reflection coefficients for the ensemble average and the periodic model.
Fig. 7 Flowchart illustrating the steps used to calculate the observed color using the periodic FDTD method. indicates the Fourier transform, and η is the intrinsic impedance of free space.
Fig. 8 Illustration for the three-dimensional representation of the observed color as a function of the direction of observation.
Fig. 9 Observed color as a function of direction for Morpho-like models. For all cases, there are eight lamellae and a base. (a) Aligned lamellae without taper, Ymax=1.0. (b) Aligned lamellae with the top four lamellae tapered, Ymax=0.30. (c) Aligned, tapered lamellae plus a cover scale (dielectric slab), Ymax=0.28. (d)–(f) Same as (a)–(c), but with offset lamellae. (d) Ymax=0.21. (e) Ymax=0.05. (f) Ymax=0.22.
Fig. 10 Far-zone scattering patterns for the Morpho model shown in Fig. 9a. (a) Pattern in the y–z plane. (b) Pattern in a perpendicular plane that passes through the peak of the beam (angle measured from θ=20°, ϕ=90°).
Fig. 11 Power reflection coefficients (a) for structures with aligned lamellae, as in Figs. 9a, 9b, 9c, and (b) for structures with offset lamellae, as in Figs. 9d, 9e, 9f.
Fig. 12 Observed color as a function of direction for Morpho-like models that include microribs. For both cases, there are eight lamellae with the top four lamellae tapered and a base: (a) aligned lamellae, Ymax=0.26; (b) offset lamellae, Ymax=0.05.
Fig. 13 Photographs of the dorsal side of a male Troides magellanus butterfly. In both views, the light source is located near the observation point. (a) From most observation directions, the observed color of the lower wings is yellow. (b) When illuminated and viewed near grazing, the color of the lower wings changes abruptly, with a change in angle, to blue.
Fig. 14 Observed color as a function of direction for Troides magellanus model with six aligned lamellae. The details for the model are the same for all cases, and only the direction of the incident light is changed: (a) θi=40°, Ymax=0.05; (b) θi=50°, Ymax=0.14; (c) θi=70°, Ymax=1.0.
Fig. 15 Power reflection coefficient for Troides magellanus model. The rapid increase in the intensity of the iridescent color is shown by comparing the reflection coefficient for θi=40° and θi=70°.
Fig. 16 Observed color as a function of direction for Ancyluris meliboeus model with six aligned lamellae and a base: (a) without microribs, Ymax=1.0; (b) with microribs, Ymax=0.48.
Fig. 17 Power reflection coefficient for Ancyluris meliboeus model. The addition of the microribs causes the peak in the reflection coefficient to shift to longer wavelengths.
1.

S. Berthier, Iridescences: The Physical Colors of Insects (Springer, 2007).

2.

T. F. Anderson and A. G. Richards, “An electron microscope study of some structural colors of insects,” J. Appl. Phys. 13, 748–758 (1942). [CrossRef]

3.

H. Ghiradella, “Structure of iridescent lepidopteran scales: variations on several themes,” Ann. Entomol. Soc. Am. 77, 637–645 (1984).

4.

H. Ghiradella, “Hairs, bristles, and scales,” in Microscopic Anatomy of Invertebrates, W. H. Frederick and L. Michael, eds. (Wiley-Liss, 1998), Vol. 11A, pp. 257–287.

5.

P. Vukusic, J. R. Sambles, and H. Ghiradella, “Optical classification of microstructure in butterfly wing-scales,” Photonics Sci. News 6, 61–68 (2000).

6.

P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wooten, “Quantified interference and diffraction in single Morpho butterfly scales,” Proc. R. Soc. London Ser. B 266, 1403–1411 (1999). [CrossRef]

7.

A. Taflove and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech, 2005).

8.

R. T. Lee and G. S. Smith, “An alternative approach for implementing periodic boundary conditions in the FDTD method using multiple unit cells,” IEEE Trans. Antennas Propag. 54, 698–705 (2006). [CrossRef]

9.

R. T. Lee, “A novel technique for incorporating periodic boundaries into the FDTD method and the application to the study of structural color of insects,” Ph.D.thesis(Georgia Institute of Technology, 2009).

10.

L. Plattner, “Optical properties of the scales of Morpho rhetenor butterflies: theoretical and experimental investigation of the backscattering of light in the visible spectrum,” J. R. Soc. Interface 1, 49–59 (2004). [CrossRef]

11.

S. Banerjee, J. B. Cole, and T. Yatagai, “Colour characterization of a Morpho butterfly wing-scale using a high accuracy nonstandard finite-difference time-domain method,” Micron 38, 97–103 (2007). [CrossRef]

12.

K. S. Yee, “Numerical solution of initial boundary value prob lems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966). [CrossRef]

13.

J. B. Schneider, “Plane waves in FDTD simulations and a nearly perfect total-field/scattered-field boundary,” IEEE Trans. Antennas Propag. 52, 3280–3287 (2004). [CrossRef]

14.

J. A. Roden and S. D. Gedney, “Convolutional PML (CPML): an efficient FDTD implementation of the CFS-PML for arbitrary media,” Microw. Opt. Technol. Lett. 27, 334–339 (2000). [CrossRef]

15.

S. Kinoshita and S. Yoshioka, “Structural colors in nature: The role of regularity and irregularity in the structure,” Chem. Phys. Chem. 6, 1442–1459 (2005). [CrossRef] [PubMed]

16.

“Standard terminology of appearance,” ASTM E 284 08 (ASTM International, 2008).

17.

R. S. Berns, Billmeyer and Saltzman’s Principles of Color Technology, 3rd ed. (Wiley, 2000).

18.

“Standard practice for computing the colors of objects by using the CIE system,” ASTM E 308 06 (ASTM International, 2006).

19.

R. A. Potyrailo, H. Ghiradella, A. Vertiatchikh, K. Dovidenko, J. R. Cournoyer, and E. Olson, “Morpho butterfly wing scales demonstrate highly selective vapour response,” Nat. Photon. 1, 123–128 (2007). [CrossRef]

20.

S. Berthier, E. Charron, and J. Boulenguez, “Morphological structure and optical properties of the wings of Morphidae,” Insect Sci. 13, 145–158 (2006). [CrossRef]

21.

S. Yoshioka and S. Kinoshita, “Wavelength-selective and anisotropic light-diffusing scale on the wing of the Morpho butterfly,” Proc. R. Soc. London Ser. B 271, 581–587 (2004). [CrossRef]

22.

B. Gralak, G. Tayeb, and S. Enoch, “Morpho butterflies wings color modeled with lamellar grating theory,” Opt. Express 9, 567–578 (2001). [CrossRef] [PubMed]

23.

C. Lawrence, P. Vukusic, and R. Sambles, “Grazing-incidence iridescence from a butterfly wing,” Appl. Opt. 41, 437–441 (2002). [CrossRef] [PubMed]

24.

J. P. Vigneron, K. Kertesz, Z. Vertesy, M. Rassart, V. Lousse, Z. Balint, and L. P. Biro, “Correlated diffraction and fluorescence in the backscattering iridescence of the male butterfly Troides magellanus (Papilionidae),” Phys. Rev. E 78, 021903 (2008). [CrossRef]

25.

P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wootton, “Limited-view iridescence in the butterfly Ancyluris meliboeus,” Proc. R. Soc. London Ser. B , 269, 7–14 (2002). [CrossRef]

26.

S. Kinoshita, S. Yoshioka, and J. Miyazaki, “Physics of structural colors,” Rep. Prog. Phys. 71, 076401 (2008). [CrossRef]

OCIS Codes
(170.1420) Medical optics and biotechnology : Biology
(330.1690) Vision, color, and visual optics : Color

ToC Category:
Vision, Color, and Visual Optics

History
Original Manuscript: March 24, 2009
Revised Manuscript: June 22, 2009
Manuscript Accepted: June 22, 2009
Published: July 14, 2009

Virtual Issues
Vol. 4, Iss. 9 Virtual Journal for Biomedical Optics
July 17, 2009 Spotlight on Optics

Citation
R. Todd Lee and Glenn S. Smith, "Detailed electromagnetic simulation for the structural color of butterfly wings," Appl. Opt. 48, 4177-4190 (2009)
http://www.opticsinfobase.org/ao/abstract.cfm?URI=ao-48-21-4177


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. S. Berthier, Iridescences: The Physical Colors of Insects (Springer, 2007).
  2. T. F. Anderson and A. G. Richards, “An electron microscope study of some structural colors of insects,” J. Appl. Phys. 13, 748-758 (1942). [CrossRef]
  3. H. Ghiradella, “Structure of iridescent lepidopteran scales: variations on several themes,” Ann. Entomol. Soc. Am. 77, 637-645 (1984).
  4. H. Ghiradella, “Hairs, bristles, and scales,” in Microscopic Anatomy of Invertebrates, W. H. Frederick and L. Michael, eds. (Wiley-Liss, 1998), Vol. 11A, pp. 257-287.
  5. P. Vukusic, J. R. Sambles, and H. Ghiradella, “Optical classification of microstructure in butterfly wing-scales,” Photonics Sci. News 6, 61-68 (2000).
  6. P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wooten, “Quantified interference and diffraction in single Morpho butterfly scales,” Proc. R. Soc. London Ser. B 266, 1403-1411(1999). [CrossRef]
  7. A. Taflove and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech, 2005).
  8. R. T. Lee and G. S. Smith, “An alternative approach for implementing periodic boundary conditions in the FDTD method using multiple unit cells,” IEEE Trans. Antennas Propag. 54, 698-705 (2006). [CrossRef]
  9. R. T. Lee, “A novel technique for incorporating periodic boundaries into the FDTD method and the application to the study of structural color of insects,” Ph.D. thesis (Georgia Institute of Technology, 2009).
  10. L. Plattner, “Optical properties of the scales of Morpho rhetenor butterflies: theoretical and experimental investigation of the backscattering of light in the visible spectrum,” J. R. Soc. Interface 1, 49-59 (2004). [CrossRef]
  11. S. Banerjee, J. B. Cole, and T. Yatagai, “Colour characterization of a Morpho butterfly wing-scale using a high accuracy nonstandard finite-difference time-domain method,” Micron 38, 97-103 (2007). [CrossRef]
  12. K. S. Yee, “Numerical solution of initial boundary value problems involving Maxwell's equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302-307 (1966). [CrossRef]
  13. J. B. Schneider, “Plane waves in FDTD simulations and a nearly perfect total-field/scattered-field boundary,” IEEE Trans. Antennas Propag. 52, 3280-3287 (2004). [CrossRef]
  14. J. A. Roden and S. D. Gedney, “Convolutional PML (CPML): an efficient FDTD implementation of the CFS-PML for arbitrary media,” Microw. Opt. Technol. Lett. 27, 334-339 (2000). [CrossRef]
  15. S. Kinoshita and S. Yoshioka, “Structural colors in nature: The role of regularity and irregularity in the structure,” Chem. Phys. Chem. 6, 1442-1459 (2005). [CrossRef] [PubMed]
  16. “Standard terminology of appearance,” ASTM E 284 08 (ASTM International, 2008).
  17. R. S. Berns, Billmeyer and Saltzman's Principles of Color Technology, 3rd ed. (Wiley, 2000).
  18. “Standard practice for computing the colors of objects by using the CIE system,” ASTM E 308 06 (ASTM International, 2006).
  19. R. A. Potyrailo, H. Ghiradella, A. Vertiatchikh, K. Dovidenko, J. R. Cournoyer, and E. Olson, “Morphobutterfly wing scales demonstrate highly selective vapour response,” Nat. Photon. 1, 123-128 (2007). [CrossRef]
  20. S. Berthier, E. Charron, and J. Boulenguez, “Morphological structure and optical properties of the wings of Morphidae,” Insect Sci. 13, 145-158 (2006). [CrossRef]
  21. S. Yoshioka and S. Kinoshita, “Wavelength-selective and anisotropic light-diffusing scale on the wing of the Morphobutterfly,” Proc. R. Soc. London Ser. B 271, 581-587 (2004). [CrossRef]
  22. B. Gralak, G. Tayeb, and S. Enoch, “Morphobutterflies wings color modeled with lamellar grating theory,” Opt. Express 9, 567-578 (2001). [CrossRef] [PubMed]
  23. C. Lawrence, P. Vukusic, and R. Sambles, “Grazing-incidence iridescence from a butterfly wing,” Appl. Opt. 41, 437-441(2002). [CrossRef] [PubMed]
  24. J. P. Vigneron, K. Kertesz, Z. Vertesy, M. Rassart, V. Lousse, Z. Balint, and L. P. Biro, “Correlated diffraction and fluorescence in the backscattering iridescence of the male butterfly Troides magellanus (Papilionidae),” Phys. Rev. E 78, 021903(2008). [CrossRef]
  25. P. Vukusic, J. R. Sambles, C. R. Lawrence, and R. J. Wootton, “Limited-view iridescence in the butterfly Ancyluris meliboeus,” Proc. R. Soc. London Ser. B , 269, 7-14 (2002). [CrossRef]
  26. S. Kinoshita, S. Yoshioka, and J. Miyazaki, “Physics of structural colors,” Rep. Prog. Phys. 71, 076401 (2008). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.


« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited