## Abstract

A novel full vector model based on the supercell lattice method is presented to analyze photonic crystal fiber (PCF). From symmetry analysis waveguide modes, we classify PCF modes into nondegenerate or degenerate pairs according to the minimum waveguide sectors and their appropriate boundary conditions. We describe how the modes of the PCF can be labeled by step-index fiber analogs, with the exception of modes that have the same symmetry as the PCF. When the doublet of the degenerate pairs both have the same symmetry as the PCF, they will be split into two nondegenerate modes.

© 2003 Optical Society of America

## 1. Introduction

Recently, intensive research has been generated by the development of photonic crystal fiber (PCF)—a new type of optical fiber that contains a central defect region in a regular lattice of air holes. According to the guidance mechanism, PCF may be divided into two general classes [1]: total-internal-reflection (TIR) PCF and photonic bandgap (PBG) PCF. Because the effective index of cladding is strongly wavelength dependent, TIR-PCF displays unusual guidance properties [2–4]: For small hole diameters, the fibers are found to be endlessly single mode; the zero dispersion point may be shifted into the visible part of the spectrum; the possibility of nearly zero ultraflattened dispersion exists; and PCF allows very strong or very weak nonlinear coefficients.

Together with the technological advancement in the fabrication of PCF, powerful theoretical tools have been developed for modeling PCF guidance properties [5–13]. These tools include the effective-index method, the plane-wave expansion method, the localized basis function method, the finite-element method, and more recently the beam propagation method and the multipole method. When higher-order modes or polarization properties are considered, the full vector approach is crucial for assessing the true behavior of electromagnetic waves in complex inhomogeneous structures such as PCF.

It is well known that PCF is often strongly multimode in the visible and near-infrared regions when the filling factor is large enough. This may lead to a number of intermodally phase-matched nonlinear processes [14]. The generation of very high-order fiber modes in the ultraviolet region has been reported [15]. It is necessary to investigate the modal properties of PCF such as degeneracy, classification, and so on. Some studies have discussed the mode properties of PCF [16–19], but the classification and degeneracy properties of higher-order modes have not been investigated adequately. In this paper, a full vector model based on our previous research is developed for modeling PCF. From the symmetry analysis of modes in waveguides, we classify the modes of the PCF into degenerate or nondegenerate mode classes according to the minimum waveguide sectors and the appropriate boundary conditions. We explain how the modes of the PCF can be labeled by its step-index fiber analogs, with the exception of modes that have the same symmetry as PCF.

## 2. Full vector model

We take the PCF to be lossless and uniform in the propagation (*z*) direction, so our main task is to investigate the transverse mode field distribution ${\overrightarrow{e}}_{t}$
(*x, y*), which can be divided into two polarization components along the *x* or the *y* axis,

where *β*_{j}
is the propagation constant of *j*th mode and *e⃗*_{t}
(*x, y*)=*e*_{x}*x̂*+*e*_{y}*ŷ*, *e⃗*_{z}
(*x, y*) are the transverse and longitudinal components of the modal electric field, respectively. Inserting Eq. (1) into the full vector wave equation, we obtain the following pair of wave equations for electric field components *e*_{x}
(*x, y*) and *e*_{y}
(*x,y*) [20],

$$\left({\nabla}_{t}^{2}-{\beta}_{j}^{2}+{k}^{2}{n}^{2}\right){e}_{y}=-\frac{\partial}{\partial y}\left({e}_{x}\frac{\partial \phantom{\rule{.2em}{0ex}}\mathrm{In}{n}^{2}}{\partial x}+{e}_{y}\frac{\partial \phantom{\rule{.2em}{0ex}}\mathrm{In}{n}^{2}}{\partial y}\right),$$

where *n*^{2}
=*n*^{2}*(x, y)* is the transverse dielectric constant profile and *k*=2π/*λ* is the wave number. To solve the coupled equations for *e*_{x}
and *e*_{y}
, the modal electric field will be expanded as

$${e}_{y}(x,y)=\sum _{a,b=0}^{F-1}{\epsilon}_{ab}^{y}{\psi}_{a}\left(x\right){\psi}_{b}\left(y\right),$$

where *F* is the number of terms retained in this expansion and *ψ*_{i}
(*s*) (*i*=*a,b, s*=*x,y*) are elements of the following orthogonal set of Hermite-Gaussian basis functions,

where *H*_{i}
(*s*/*ω*) is the *i*th-order Hermite polynomial and *ω* is the characteristic width of the basis set.

We use the supercell lattice method [21] for describing the transverse dielectric structure of PCF. To describe the dielectric structure of the supercell lattice, two 2-D virtually perfect photonic crystals (PC1 and PC2) are introduced, which can be added together to reconstruct the dielectric structure of the supercell lattice. The dielectric constant of the PCF can be expressed as the sum of the decompositions of both PC1 and PC2:

where *P*
_{1} and *P*
_{2} are the number of decomposition terms of PC1 and PC2, respectively; *l*
_{1x}, *l*
_{1y}, *l*
_{2x}, and *l*
_{2y} are the characteristic lengths along the *x* or the *y* axis of PC1 and PC2; and *P*
_{1ab}, *P*
_{2ab}, ${P}_{\mathrm{l}\mathit{\text{ab}}}^{1\mathrm{n}}$, and ${P}_{2\mathit{\text{ab}}}^{1\mathrm{n}}$ are the decomposition coefficients that can be obtained through the Fourier transform [21, 22] of the periodic dielectric structures for both PC1 and PC2.

All decomposition coefficients in Eq. (5) can be evaluated analytically. The supercell overlapping method can be used to analyze many kinds of PCF with different lattice structures. Figure 1 demonstrates the simulation results of the dielectric structure of TIR-PCF. It is obvious that all these simulation results can reflect the dielectric structure accurately, so it is shown that the supercell overlapping method is accurate and reliable.

Substituting the decompositions Eq. (3) into the wave Eqs. (2) by use of the orthonormality of the Hermite-Gaussian basis functions, the eigenvalue equation will be obtained and shown as Eq. (6),

where *I*
^{(1)}, *I*
^{(2)} and *I*
^{(3)}, *I*
^{(4)} are overlapping integrals of the modal functions; they are defined as

$${I}_{\mathit{abcd}}^{\left(2\right)}={\iint}_{-\infty}^{+\infty}{n}^{2}{\psi}_{a}\left(x\right){\psi}_{b}\left(y\right){\psi}_{c}\left(x\right){\psi}_{d}\left(y\right)\mathit{dx}\mathit{dy},$$

$${I}_{\mathit{abcd}}^{\left(3\right)x}={\iint}_{-\infty}^{+\infty}{\psi}_{a}\left(x\right){\psi}_{b}\left(y\right)\frac{\partial}{\partial x}\left[{\psi}_{c}\left(x\right){\psi}_{d}\left(y\right)\frac{\partial \phantom{\rule{.2em}{0ex}}\mathrm{In}{n}^{2}}{\partial x}\right]\mathit{dx}\mathit{dy},$$

$${I}_{\mathit{abcd}}^{\left(3\right)y}={\iint}_{-\infty}^{+\infty}{\psi}_{a}\left(x\right){\psi}_{b}\left(y\right)\frac{\partial}{\partial y}\left[{\psi}_{c}\left(x\right){\psi}_{d}\left(y\right)\frac{\partial \phantom{\rule{.2em}{0ex}}\mathrm{In}{n}^{2}}{\partial y}\right]\mathit{dx}\mathit{dy},$$

$${I}_{\mathit{abcd}}^{\left(4\right)x}={\iint}_{-\infty}^{+\infty}{\psi}_{a}\left(x\right){\psi}_{b}\left(y\right)\frac{\partial}{\partial x}\left[{\psi}_{c}\left(x\right){\psi}_{d}\left(y\right)\frac{\partial \phantom{\rule{.2em}{0ex}}\mathrm{In}{n}^{2}}{\partial y}\right]\mathit{dx}\mathit{dy},$$

$${I}_{\mathit{abcd}}^{\left(4\right)y}={\iint}_{-\infty}^{+\infty}{\psi}_{a}\left(x\right){\psi}_{b}\left(y\right)\frac{\partial}{\partial y}\left[{\psi}_{c}\left(x\right){\psi}_{d}\left(y\right)\frac{\partial \phantom{\rule{.2em}{0ex}}\mathrm{In}{n}^{2}}{\partial x}\right]\mathit{dx}\mathit{dy}.$$

Substituting decomposition Eq. (5) into Eq. (7), we can calculate the overlap integrals by using the orthonormality of the Hermite-Gaussian basis functions and the standard integration-by-parts technique combined with some of the definite integrals available in collections [23] along with some identities [24]. The fact that all these integrals can be calculated analytically is a significant advantage: The numerical precision and efficiency can be improved greatly, especially when a large number of terms are needed in the expansions. Once evaluated, the expressions for these overlapping integrals become quite complicated; they are not given here, but the analogous expressions can be found in Refs. [21] and [22].

When we solve Eq. (6) at a particular wavelength *λ*, the modes and the corresponding *β*_{j}
of the PCF at *λ* can be calculated. Because * L* is a (2×

*F*

^{2})×(2×

*F*

^{2}) matrix, the eigenvalue problem can be solved to produce 2×

*F*

^{2}eigenvalues and their corresponding eigenvectors. The guided mode solutions can be distinguished from the radiation modes by extraction of the eigenvalues, which fall within the following range allowed by the structure:

*n*

_{eff}<

*β*/

*k*<

*n*

_{si}, where

*n*

_{eff}is the effective index [2,21] of the photonic crystal cladding at the wavelength

*λ*. Only a few of the eigenvalue-eigenvector pairs are physically relevant quantities corresponding to the guided modes of the fiber. The other eigenvalues are related to radiation modes of the fiber. The transverse modal electric field can be obtained by means of substituting the eigenvector into decomposition Eq. (5).

## 3. Modes in PCF

In scalar approximation, the characteristics of polarization in the PCF are hidden. With the full vector approach developed in Section 2, we can predict the vector modal behavior in PCF. In this section, we discuss the classification of modes in the PCF according to the analysis of symmetry. Compared with step-index fiber, the modes of the PCF are labeled by its analogs.

#### 3.1 Analysis of symmetry

From the theory of group representations, we know that the symmetry of a waveguide controls several important characteristics of the modes of the waveguide [25]. Determining the symmetry type of a particular waveguide enables one to classify the possible modes in mode classes, predict the mode degeneracies between mode classes, and determine the azimuthal symmetries of modal electromagnetic fields in each mode class. Further, from the azimuthal symmetry the modal electromagnetic fields of a mode class, one can specify a minimum sector of waveguide cross section, together with associated boundary conditions for this sector, which is necessary and sufficient for completely determining all the modes of that mode class.

If a uniform waveguide has *n*-fold rotation symmetry and also possesses precisely *n* reflection planes spaced azimuthally by π/*n* radians, it is said to have C_{nν} symmetry. Because the cross section of the triangular lattice PCF studied here has sixfold rotation symmetry in addition to mirror symmetry, the point group is C_{6ν}. The C_{6ν} symmetry of PCF suggests the following: All the modes can be divided into eight classes; these classes either exhibit the full waveguide symmetry and are nondegenerate, or the occur in degenerate pairs that support this symmetry only in combination [12,18,19].

Figure 2 shows the minimum sector for waveguides with C_{6ν} symmetry. Mode classes *p*=1,2,7,8 are nondegenerate, and they exhibit the full symmetry of the PCF structure. *p*=3,4 and *p*=5,6 are degenerate pairs, and they exhibit the full symmetry of structure in combination. The difference between the modal classes is different minimum sector combining with boundary conditions applied to the tangential component of electric or magnetic fields at the edge of the minimum sector. These boundary conditions are either short circuit (*e*_{z}
=0) or open circuit (*h*_{z}
=0), or they are combinations of these (see Fig. 2). The use of the minimum waveguide sector for analyzing different mode classes leads to an important advantage: For mode classes *p*=1,2,7,8, with proper boundary conditions (see Fig. 2), only the π/6 sector is needed for calculating the modal field, and only the π/2 sector is needed for *p*=3,4 and *p*=5,6. This approach is used in the multipole method for modeling PCF [12,13]. In this paper, we will check the mode classification and degeneracy rather than take advantage of the symmetry analysis for calculating modal field.

It is well known that the transverse electric and magnetic fields in a uniform waveguide can be expressed in terms of the longitudinal components. Equation (8) shows the transverse electric field as an example [20]:

McIsaac [25] expressed longitudinal electric- and magnetic-field components in terms of a Fourier series in the azimuthal angle *θ* for uniform waveguides with C_{6ν} symmetry. Substituting these expressions into Eq. (8), we find that the short-circuit boundary conditions (*e*_{z}
=0, solid lines in Fig. 2) are equal only to azimuthal components existing in transverse electric-field expressions, namely, *e⃑*_{t}
=*$\widehat{\phi}$
e _{φ}*. In a similar manner, open-circuit boundary conditions (

*h*

_{z}=0, dashed lines in Fig. 2) are equal only to radial components existing in the transverse electric-field expression, namely,

*e⃑*

_{t}=

*r̂e*

_{r}.

From this analysis, we can classify the modes of the PCF by the behavior of transverse electric fields at the boundary of each minimum sector in Fig. 2.

The higher-order modes will exist only if the relative hole size *d*/Λ is greater than a certain value [3,5]. We apply our approach to TIR-PCF with structural parameters *Λ*=2.3µm, *d*/*Λ*=0.8, and wavelength 633nm. From Eq. (6), we extract 24 eigenvalues *β*_{j}
/*k* below *n*_{si}
and plot them in Fig. 3. Figure 3 shows the mode indices of different modes in TIR-PCF. We find that only the first 14 modes are bound, and the other 10 are radiation modes. We list the first 14 modes with mode index (*n*_{eff}
), mode class (*p*), degeneracy, computation error (*Δn*), and label in Table 1. We discuss these modes in the following sections.

#### 3.2 Fundamental mode

Extracting the first two eigenvalue-eigenvector pairs, we obtain the fundamental polarization doublet modes in PCF. Figure 4 shows the total modal intensity distribution and 2-D electric vector distributions of the polarization doublet modes. The refractive-index profile is also superposed in Fig. 4. According to the 2-D electric vector distributions, we find at θ=0 that only azimuthal or radial components of the electric field exist for Figs. 4(b) and 4(c); the case is reversed at θ=π/2. So the components belong to mode classes *p*=4 and 3, respectively. These two modes must be fully degenerate theoretically. Figure 4(a) shows a combined intensity distribution of the doublet modes.

The electric field distributions looks exactly the same as the HE_{11} mode in step-index fiber, and we can label them by the degenerate pair of the HE_{11} mode [17]. We find that there is a slight difference between the electric field distributions of these two modes. According to the symmetry analysis, even though the *x* and *y* directions are not equivalent for the triangular lattice PCF, their mode index must be the same theoretically. As a result, the observation of the birefringence must be a result of the error produced by the numerical calculations. When an ideal PCF is examined, the modal birefringence Δ*n*=*n*_{x}*-n*_{y}
can be used to scale the accuracy of the algorithm [21]. The less |Δ*n*| is, the more accurate the numerical method. The modal birefringence Δ*n*=*n*_{x}*-n*_{y}
is of the order of 10^{-7}, which scales the accuracy of the algorithm.

#### 3.3 First higher-order modes

In Fig. 3, the mode indices from the third to the sixth mode are almost same; they belong to the first higher-order multiplet. We show electric vector distributions for these modes in Fig. 5. The four modes of the first higher-order multiplet are almost degenerate, and the mode indices are 1.43674688, 1.43637588, 1.43635241, and 1.43622157, respectively. The electric vector distributions in Figs. 5(a) and 5(d) resemble the TE_{01} and TM_{01} modes in step-index fiber, and we can label them by the TE_{01} and TM_{01} modes. From the mode classification in Fig. 2, only the azimuthal component of the electric vector in Fig. 5(a) exists at angle 0 and π/6, and only the radial component of the electric vector in Fig. 5(d) exists at angle 0 and π/6; they belong to classes *p*=2 and *p*=1, respectively. Moreover, these two modes are related in a sixfold rotationally invariant fashion to the triangular lattice PCF, so they are nondegenerate. The mode-index difference between TE_{01} and TM_{01} modes is 5.253053×10^{-4}, it cannot be considered to be a computational error. The other two modes [Figs. 5(b) and 5(c)] belong to classes *p*=5 and *p*=6 because the electric vector distributions satisfy corresponding boundary conditions at θ=0 and π/2. According to the electric vector distributions, they can be considered to be a degenerate pair of the HE_{21} mode with computation error (mode-index difference) of 2.346547×10^{-5}. In Figs. 5(b) and 5(c), the electric vector field distributions do not exhibit the symmetry of the triangular lattice PCF, so they must provide degenerate pair support of the symmetry of PCF in combination. Figure 6 shows the combination of modes 4 and 5; the intensity distribution reflects the C_{6ν} symmetry of the PCF. The intensity distributions of the other two modes (3 and 6) are not shown because they are similar to the HE_{21} mode.

#### 3.4 Other modes

Figure 7 shows electric vector distributions of modes 7–14 in Table 1. According to the electric vector distributions, we can divide them into different classes according to the symmetry analysis. These classifications are listed in Table 1, and we label them by their step-index fiber analogs. There is an interesting phenomenon here: for modes 7 and 10, their electric vector distributions look like the HE_{31} mode of the step-index fiber, but they are nondegenerate from the symmetry analysis. The electric vector distributions for modes 7 and 10 can exhibit full waveguide symmetry, and the degeneracy of the HE_{31} mode is split. This conclusion is also supported by their mode-index difference 2.307674×10^{-3}, which is much larger than the computation error, such as the mode-index difference 4.396536×10^{-5} between modes 8 and 9. In Fig. 8 we present the intensity distributions of the nondegenerate mode HE_{311} and the degenerate modes EH_{11}, HE_{12}, and EH_{21}. The intensity distributions of mode HE_{312} are very similar to those of mode HE_{311} except for the amplitude, which is not shown here. These intensity distributions exhibit C_{6ν} symmetry related to the structure of PCF.

From the analysis of the PCF modes, it can be deduced that the degeneracy of the higher-order mode analogs of the step-index fiber with C_{6ν} symmetry will be split into triangular lattice PCF, for example, HE_{31} mode and HE_{3s,t} EH_{3s, t} modes, where *s* and *t* are integers. In fact, when relative hole size *d*/*Λ* is 0.9 (the other parameters fixed), modes 21 and 22 are bound, and they can be labeled by EH_{311} and EH_{312} with a mode-index difference of the order of 10^{-4}, which is 1 order of magnitude larger than the computation error. Figure 9 is a modal intensity distribution of mode 21 and 2-D electric vector distributions of modes 21 and 22, which belong to classes 8 and 7, respectively. Our conclusion is different from that of Ref. [19], in which all HE-like modes are considered to be degenerate.

## 4. Discussion

Symmetry analysis can provide information about the general characteristics of the modes of waveguides. We can classify all the modes into mode classes on the basis of the modal field symmetry, and we can determine the possible mode degeneracies between mode classes and specify the minimum sector of the waveguide cross section with its appropriate boundary conditions, which can completely determine the modes of a given mode class.

The modes of a waveguide with C_{nν} symmetry either exhibit the full waveguide symmetry (C_{nν}) and are nondegenerate, or they occur in degenerate pairs that support this symmetry only in combination. When PCF is examined, the doublet of the degenerate pairs in which both have the same symmetry as PCF will be split into two nondegenerate modes. For instance, the degenerate mode HE_{ns/2,t} or EH_{ns/2,t} (*s* and *t* are integers; *n* is even) of waveguides with C_{nν} symmetry will be split into two nondegenerate modes.

For the triangular lattice PCF, as discussed above, the degeneracy will be split, because each mode of the mode pairs with C_{6ν} symmetry can reflect the full waveguide symmetry.

For the square lattice PCF with C_{4ν} symmetry, Fig. 10 shows the minimum sectors for waveguides with C_{4ν} symmetry. Mode classes *p*=1,2,5,6 are nondegenerate; they exhibit full the symmetry of the PCF structure. Modes *p*=3,4 are degenerate pairs, and they exhibit the full symmetry of the structure in combination. Table 2 lists mode index, degeneracy, mode class, computation error, and label of the first eight bound modes of the square lattice PCF with parameters *Λ*=2.3µm, *d*/*Λ*=0.8 and wavelength 633nm.

The HE_{21} mode is split into HE_{211} and HE_{212} mode (see Table 2) because the modal field exhibits the same symmetry (C_{4ν}) as that of the square lattice PCF. Figure 11 shows the electric vector distributions of modes 3–6 (HE_{211}, TM_{01}, TE_{01}, and HE_{212}), and each of these four modes reflects C_{4ν} symmetry and is nondegenerate. Figure 12 gives the intensity distributions of mode 3 (HE_{211}) and mode 5 (TE_{01}), and the other two modes have intensity distributions similar to those of mode 3, which are not shown here. In the same way, it is deduced that the degeneracy of mode with C_{4ν} symmetry will be split, for example, into HE_{2s,t}, EH_{2s,t} modes.

## 5. Conclusion

A full vector method for modeling PCF has been presented. On the basis of the novel supercell lattice method, the transverse-index profile of PCF is expanded in cosine functions. The propagation constant and the mode field distribution of the PCF can be calculated by means of recasting the Maxwell equations into an eigenvalue system. With this new method, the classification and degeneracy of modes in PCF are discussed on the basis of symmetry analysis. The modes of a waveguide with C_{nν} symmetry either exhibit the full waveguide symmetry (C_{nν}) and are nondegenerate or occur in degenerate pairs that support this symmetry only in combination. We group the modes of the PCF into degenerate or nondegenerate mode classes according to the minimum waveguide sectors and appropriate boundary conditions. We show that the modes of the PCF can be labeled by the step-index fiber analogs. When mode pairs have the same symmetry as PCF, they will be split into two nondegenerate modes. It is deduced that the degenerate mode HE_{ns/2,t} or EH_{ns/2,t} of PCF with C_{nν} symmetry will be split into two nondegenerate modes.

## References and links

**1. **J. C. Night and P. St. Russell, “New way to guide light,” Science **296**, 276–277 (2002). [CrossRef]

**2. **T. A. Birks, J. C. Knight, and P. St. J. Russell, “Endlessly single-mode photonic crystal fiber,” Opt. Lett. **22**, 961–963 (1997). [CrossRef] [PubMed]

**3. **A. Ferrando and E. Silvestre, “Nearly zero ultraflattened dispersion in photonic crystal fibers,” Opt. Lett. **25**, 790–792 (2000). [CrossRef]

**4. **A. L. Gaeta, “Nonlinear propagation and continuum generation in microstructured optical fibers,” Opt. Lett. **27**, 924–926 (2002). [CrossRef]

**5. **T. A. Birks, D. Mogilevtsev, J. C. Knight, P. St. J. Russell, J. Broeng, P. J. Roberts, J. A. West, D. C. Allan, and J. C. Fajardo, “The analogy between photonic crystal fibres and step index fibres,” in *Optical Fiber Communication Conference* (Optical Society of America, Washington, D.C., 1998), FG4, pp. 114–116.

**6. **J. D. Joannopoulos, R. D. Meade, and J. N. Winn, *Photonic Crystals: Molding the Flow of Light* (Princeton U. Press, New York, 1995).

**7. **S. Guo and S. Albin, “Simple plane wave implementation for photonic crystal calculations,” Opt. Express **11**, 167–175 (2003), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-11-2-167 [CrossRef] [PubMed]

**8. **T. M. Monro, D. J. Richardson, N. G. R. Broderick, and P. J. Bennett, “Modeling large air fraction holey optical fibers,” J. Lightwave Technol. **18**, 50–56 (2000). [CrossRef]

**9. **D. Mogilevtsev, T. A. Birks, and P. St. Russell, “Localized function method for modeling defect mode in 2-d photonic crystal,” J. Lightwave Technol. **17**, 2078–2081 (1999). [CrossRef]

**10. **M. Koshiba, “Full vector analysis of photonic crystal fibers using the finite elelment method,” IEICE Electron. E85-C , **4**, 881–888 (2002).

**11. **F. Fogli, L. Saccomandi, P. Bassi, G. Bellanca, and S. Trillo, “Full vectorial BPM modeling of index-guiding photonic crystal fibers and couplers,” Opt. Express **10**, 54–59 (2002), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-10-1-54 [CrossRef] [PubMed]

**12. **T. P. White, B. T. Kuhlmey, R. C. McPhedran, D. Maystre, G. Renversez, C. Martijn de Sterke, and L. C. Botten, “Multipole method for microstructured optical fibers. I. Formulation,” J. Opt. Soc. Am. B **19**, 2322–2330 (2002). [CrossRef]

**13. **B. T. Kuhlmey, T. P. White, G. Renversez, D. Maystre, L. C. Botten, C. Martijn de Sterke, and R. C. McPhedran, “Multipole method for microstructured optical fibers. II. Implementation and results,” J. Opt. Soc. Am. B **19**2331–2340 (2002). [CrossRef]

**14. **F. G. Omenetto, A. J. Taylor, M. D. Moores, J. Arriaga, J. C. Knight, W. J. Wadsworth, and P. St. J. Russell, “Simultaneous generation of spectrally distinct third harmonics in a photonic crystal fiber,” Opt. Lett. **26**, 1158 (2001). [CrossRef]

**15. **A. Efimov, A. J. Taylor, F. G. Omenetto, J. C. Knight, W. J. Wadsworth, and P. St. J. Russell, “Nonlinear generation of very high-order UV modes in microstructured fibers,” Opt. Express **11**, 910–918 (2003), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-11-8-910 [CrossRef] [PubMed]

**16. **F. Brechet, J. Marcou, D. Pagnoux, and P. Roy, “Complete analysis of the characteristics of propagation into photonic crystal fibers, by the finite element method,” Opt. Fiber Technol. **6**, 181–191 (2000). [CrossRef]

**17. **A. Ferrando, E. Silvestre, J. J. Miret, P. Andrés, and M. V. Andrés, “Vector description of higher-order modes in photonic crystal fibers,” J. Opt. Soc. Am. B **17**, 1333–1340 (2000). [CrossRef]

**18. **M. J. Steel, T. P. White, C. M. Sterke, R. C. McPhedran, and L. C. Botten, “Symmetry and degeneracy in microstructured optical fibers,” Opt. Lett. **26**, 488–490 (2001). [CrossRef]

**19. **M. J. Steel, T. P. White, C. M. Sterke, R. C. McPhedran, and L. C. Botten, “Symmetry and birefringence in microstructured optical fibers,” *Optical Fiber Communication Conference 2001* (Optical Society of America, Washington, D.C., 2001), WDD88-1-WDD88-3.

**20. **A. W. Snyder, *Optical Waveguide Theory* (Chapman and Hall, New York, 1983).

**21. **W. Zhi, R. Guobin, L. Shuqin, and J. Shuisheng, “Novel supercell lattice method for the photonic crystal fibers,” Opt. Express **11**, 980–991 (2003), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-11-9-980. [CrossRef] [PubMed]

**22. **W. Zhi, R. Guobin, L. Shuqin, and J. Shuisheng are preparing a manuscript to be called “A novel supercell overlapping method for different photonic crystal fibers.”

**23. **I. S. Gradshtein and I. M. Ryzhik, *Tables of Integrals, Series and Products* (Academic, New York, 1994).

**24. **I. Kimel and L. R. Elias, “Relations between Hermite and Laguerre Gaussian modes,” IEEE J. Quantum Electron. **29**, 2562–2567 (1993). [CrossRef]

**25. **P. R. McIsaac, “Symmetry-induced modal characteristics of uniform waveguides. I. Summary of results,” IEEE Trans. Microwave Theory Tech. **MTT-23**, 421–429 (1975). [CrossRef]