## Abstract

We present a useful framework based on the coupled-wave theory, through which we can survey the resonant modes of TM polarization in 2D photonic-crystal lasers and understand their properties in detail. Through numerical calculations, we clarify their threshold gains, deviations from the Bragg frequency and field distributions. We find that the lasing mode can be selected by manipulating the hole-filling factor or the boundary reflection.

©2007 Optical Society of America

## 1. Introduction

Two dimensional (2D) photonic crystal (PC) lasers [1–6] have been attracting much attention because of their important characteristics such as single mode operation in a large cavity area and a diffraction-limited, circular-shaped output beam. The former can be achieved due to the 2D distributed-feedback (DFB) mechanism, while the latter is attributed to the surface emission from an isotropic lasing area. Furthermore, the PC lasers can tailor the shape of the output beam by controlling the shape of the unit cell or the periodicity in the PC structure [4]. Most of these experimental data agreed well with a theoretical analysis using the 2D plane wave expansion method (PWEM) [3, 7, 8] and the finite difference time domain (FDTD) method [7, 9]. However these theoretical analyses have some limitations. Two dimensional PWEM only applies to infinite structures and the FDTD method consumes enormous computer resources to analyze the actual finite structure. Hence the previous FDTD method approximated the model with much fewer periods than the actual structure.

A coupled wave treatment developed by Kogelnik and Shank (KS) [10] can circumvent such limitations. In the very limited work to date on the coupled wave analysis of 2D PC lasers [6], the discussion has been confined to a very special case, in which lasing oscillation on symmetric modes can be realized. However, novel shapes of the output beam are attributed to anti-symmetric modes [4]. Hence, it is also important to investigate the detailed properties of various resonant modes including anti-symmetric modes. Specifically, the mode spectrum (plot of threshold gain against frequency deviation from the Bragg condition), the mode pattern (intensity distribution throughout the PC area) and the effect of external reflection are very important aspects of PC laser performance. In this paper, we extend the coupled wave theory to a 2D PC structure keeping the simple form of KS’s work, and develop a comprehensive understanding of the 2D resonant modes for square lattice structure with TM polarization. It is important to understand the resonant modes with TM polarization as well, since the promising light sources for THz region, e.g. quantum-cascade lasers, have their gain in TM polarization [11]. In addition, the surface emission from a particular anti-symmetric mode of TM polarization forms a radially polarized doughnut-shaped beam [4], which can be focused to a size much less than the wavelength [12]. In section two, we derive coupled wave equations from a four wave model. In section three, numerical results for the resonant modes and the dependence on the coupling constant, κ, are presented. In section four, effects of the boundary reflections are investigated.

## 2. Coupled-wave model and equations

The 2D DFB structure investigated here consists of a square lattice with circular holes in the *x* and *y* directions, as shown in Fig. 1. The structure is assumed to be uniform in the *z* direction.

The scalar wave equation for the electric field *E _{z}* in the TM mode is written as [13],

Where [14]

*β*=*n*
_{av}
*ω*/*c* where *n*
_{av} is the averaged refractive index, *α* is an averaged gain in the medium, **G**=(m*β*
_{0}, n*β*
_{0}) is the reciprocal lattice vector, m and n are arbitrary integers, *β*
_{0}=2*π*/*a* where *a* is the lattice constant and *κ*(**G**) is the coupling constant. In the derivation of Eq. (2), we have set *α*, *α*
_{G}≪*β*
_{0} and *n*
_{G}≪*n*
_{0}. In Eq. (2), the periodic variation in the refractive index is included as a small perturbation and appears in the third term through the coupling constant *κ*(**G**) of the form,

where *n*(**G**) is the Fourier coefficient of the periodic refractive-index modulation and *λ* is the Bragg wavelength given by *λ*=*a n*
_{av}. In principle, a periodic perturbation generates an infinite set of diffraction orders. However, in the vicinity of the Bragg frequency, only a few diffraction orders have significant contributions. We have to carefully choose the appropriate diffraction orders for the Bragg frequency under investigation. Note that in this paper we consider the Bragg frequency corresponding to the Γ point in the photonic bandstructure [8], in which a second order grating induces 2D optical coupling and surface emission.

In an infinite periodic structure, the electric field is given by the Bloch mode [13],

where *h*
_{G} is the amplitude of each plane wave and **k** is a wave vector in the first Brillouin zone which becomes zero at the Γ point. For a finite structure, the amplitude of each plane wave is not a constant, so *h*
_{G} becomes a function of space. In the present work, we include four basic waves propagating in the Γ-X direction with (|**k**+**G**|=*β*
_{0}). The contributions of other waves of higher order with (|**k**+**G**|≥ √2 *β*
_{0}) in the Bloch mode are considered to be negligible. Therefore, four waves are considered in the 2D coupled-wave model. The electric field can be rewritten as,

where *R _{x}*,

*S*,

_{x}*R*and

_{y}*S*are the amplitudes of the four basic waves propagating in the +

_{y}*x*, -

*x*, +

*y*, -

*y*directions, respectively. These amplitudes correspond to

*h*

_{G}in Eq. (4). We should note here that, in contrast to TM polarization, a model for TE polarization should include higher order waves because of the different coupling mechanism [15]. All of the four amplitudes are assumed to vary slowly so that their second spatial derivatives in Eq. (1) can be neglected. We insert Eq. (2) and Eq. (5) into Eq. (1), compare terms with the same exponentials and obtain a set of coupled wave equations of the form,

In order to include surface-emission loss, we have added a surface-coupling coefficient κ_{0} [6, 16]. The frequency deviation is defined as *δ*= (*β*
^{2}-*β*
_{0}
^{2})/2*β* ≈ *β*-*β*
_{0}. As the frequency deviations are assumed to be small, we have set *β*/*β*
_{0} ≈ 1 in the above derivation. We have to note here that *δ* or *β* is an eigenvalue for each resonant mode. In Eq. (6), the coupling constants are rewritten as,

Note that the coupling constants κ_{2} and κ_{3} indicate the coupling in orthogonal and backward direction, respectively. The coupling constant κ_{1}, which induces the surface emission, is assumed to be proportional to *κ*
_{0}. The values of these constants for circular holes will be derived in the next section.

Unlike the coupled-wave equations in one dimension, the general solutions for Eq. (6) are not easily found. In order to numerically solve Eq. (6), we have used the finite-difference method to obtain an eigenvalue problem. The eigenvalue (α-*i*δ) can be obtained by solving the general eigenvalue problem using Linear Algebra PACKage (LAPACK) for the given coupling constants, κ_{2} and κ_{3}.

## 3. Numerical results

The numerical results are obtained by solving the complex simultaneous equations for the general eigenvalues (α-*i*δ). We present the data for the threshold gains, frequencies (wavenumber deviation from the Bragg condition) and intensity patterns of the resonant modes of the 2D PC structure. The calculating area is a square of the length L. The facet reflectivity is neglected in this section, but the next section describes its effect. The filling factor (area of a hole / area of the unit cell) of the circular holes is set to about 16%, where the ratio of the two coupling constants κ_{2}/κ_{3} was 2. We only consider the case with index coupling.

#### 3.1. Threshold gain and frequency

Figure 2 shows the numerical results for the resonances in the case of κ_{2}L=-4 and κ_{3}L=-2 (the reasoning will be given in section 4). The surface-coupling coefficient was set to be κ_{0}L=0.09 (which was estimated using the value in Ref. [6]). Here we have plotted the threshold gain (αL) as a function of the wavenumber deviation from a Bragg condition (δL) in Fig. 2(a). Each point represents an individual mode and is categorized by the order of the longitudinal mode pattern (the intensity profile along the direction of the wave propagation) as shown in Fig. 2(b).

If we set κ_{2}L to zero, removing the two dimensional coupling, the resonant modes become degenerate and identical to the 1D DFB results of KS’s analysis. The order (1st, 2nd, 3rd) corresponds to their N number in Fig. 5 of Ref. 10. The differences within a given order will be explained in the next section in terms of the transverse mode pattern.

#### 3.2. Mode pattern

Figure 2(b) shows the intensity distribution, |*R _{x}*|

^{2}+ |

*S*|

_{x}^{2}, along their propagating direction (

*x*), e.g. the longitudinal mode pattern for 1st, 2nd and 3rd modes. The intensity distribution along the axis (

*y*) perpendicular to the propagating direction is the transverse mode pattern, which distinguishes the individual modes within an order. Figure 3(a) is a plot of the 1st order modes extracted from Fig. 2(a) with the point labeled. We have named the basic modes A, B, E, O and O’, while labeling the higher order modes 1-4, 1’-4’. Figure 3(b) shows the amplitude sum,

*R*+

_{x}*S*, for modes A, B, E, O, O’, 1-4 and 1’-4’. Even though their longitudinal order is the same, the basic modes A, B, E, O and O’ have the lowest transverse order, while the modes 1-4 and 1’-4’ have increasingly higher transverse orders. Note here that there are certain amount of modes with the intermediate profiles between modes 3 (3’) and 4 (4’).

_{x}The amplitude sum of the two waves propagating in the *y* direction, *R _{y}* +

*S*, also shows that the basic modes A, B, E, O, O’ have the lowest transverse mode patterns. Therefore, the intensity distribution of all the waves, |

_{y}*R*|

_{x}^{2}+|

*S*|

_{x}^{2}+ |

*R*|

_{y}^{2}+|

*S*|

_{y}^{2}, for the basic modes show the lowest mode profile as shown in Fig. 4(a). Modes A and B have the single-lobed pattern and mode E (degenerate) has a saddle shaped pattern. These patterns are due to the symmetric nature of the longitudinal and transverse mode profile. On the other hand, modes O and O’ have a basin-shaped pattern with zero intensity at the center due to the anti-symmetric nature of the mode profile. The electric field

*E*given by Eq. (5) around the center of the cavity is shown in Fig. 4(b), which indicates modes A, B and E have the A

_{z}_{1}, B

_{1}and E symmetry [9], respectively. Hence modes A, B and E are considered to correspond with the resonant modes at the band edge. The surface emissions from modes A and B have an anti-symmetric nature in the near field, while mode E has a symmetric nature in the near field [6]. We have to note here that mode A produces a novel beam profile: radially-polarized, doughnut-shaped beam. This beam profile is very important for light sources with super-resolution [4], which can be focused to a size much less than the wavelength [12]. According to the intensity profile of the mode pattern, selection between mode A and mode O can be achieved by different excitation conditions: Mode A can oscillate when the excitation is made for the central part of the device, while mode O can oscillate when the excitation is made for the edge. Such mode selection is also found to be achievable by controlling the facet reflectivity, which will be discussed in section 4.

#### 3.3. Dependence on the coupling constants

The product of the coupling constant and the cavity length, κL, is the important parameter for the 1D DFB lasers. For example, in order to avoid spatial hole-burning, the intensity distribution in a cavity should be uniform, which requires κL to be around 1–4. Here we have investigated the dependence of the threshold gain and mode profile on the coupling constants.

The coupling constants are functions of the index, the gain variation, and the hole-filling factor. As shown in Eq. (3), the coupling constants are determined by the Fourier coefficient of the periodic modulation for the refractive index, *n*(**G**), and the gain, *α*(**G**). For the square lattice of circular holes, these coefficients are written as [13],

where, *n _{a}* and

*n*are the respective refractive indices for the hole and the surrounding material,,

_{b}*α*and

_{a}*α*are the gain for these two elements and

_{b}*J*

_{1}(

*x*) is the Bessel function of the first order. Substituting Eq. (8) and (9) in to Eq. (3), we obtained the coupling constant

*κ*(

**G**). Figure 5 shows the product of the coupling constants and the cavity length for the index coupling, where

*α*(

**G**)=0 in the case of

*n*<

_{a}*n*. Over

_{b}*f*≈0.79, the holes are connected to one another. At

*f*≈ 0.16, the ratio of the two coupling constants is estimated to be κ

_{2}L/κ

_{3}L=2.

We set the ratio of the two coupling constants to κ_{2}L/κ_{3}L=2 and changed their magnitude. Figure 6(a) shows the threshold gain and the wavenumber deviation of the lowest gain mode A. As the coupling constant increases, the threshold gain decreases and the wavenumber deviation increases. Figure 6(b) shows the intensity distribution of the four waves, |*R _{x}*|

^{2}+|

*S*|

_{x}^{2}+|

*R*|

_{y}^{2}+|

*S*|

_{y}^{2}, for mode A along the

*x*axis at the center of the cavity. The data were normalized by the value at the cavity edge (

*x*= ±L/2). Unlike the 1D DFB case, we could not see any obvious changes in the patterns. The intensity peaks at the center and decays toward the edges even for the small value of κ

_{3}L=0.8. This may be because κ

_{2}L=2 still conserves the feedback mechanism in which the waves couple to one another to keep the total energy inside the cavity. If the value of κ

_{3}L increases, the intensity at the center increases since the energy is more strongly confined to the center area.

Next we investigate the dependence on the hole-filling factor. Figure 7 plots the threshold gain of the modes A, B, E as functions of the hole-filling factor. The solid curves are calculated using the spline interpolation.

Among all the resonant modes, mode A has the lowest threshold gain for *f* < 0.3 and *f* > 0.6, while mode B has the lowest for 0.3 < *f* < 0.6. This result shows that the hole-filling factor controls the mode selection scheme.

## 4. Boundary reflection

Until now, we have neglected the effect of the boundary reflector. However, it is natural to have the refractive index change at the boundary of the 2D PC area. Hence the boundary reflector should exist and have some effect on the resonant mode. As was already shown for 1D DFB lasers, threshold gain was determined by the reflectivity and the phase at the end of the cavity [17, 18]. Here we investigate the external reflection at the boundary of the 2D PC cavity.

The reflection at the facet is written as,

where ρ is the reflectivity and *ϕ* is the total phase. These two parameters are the functions of position, *x* and *y*. The total phase *ϕ* is dependent on two things: the edge’s position relative to the peaks and troughs of the grating, and whether the PC edge connects to an area of higher or lower refractive index. Figure 8 shows the phase *ϕ* for a typical cleaved facet in the case of the second order grating.

The new boundary conditions at the facet (*x*=±*L*/2 or *y*=±*L*/2) for the four waves become,

where ρ_{x±} (ρ_{y±}) and *ϕ*
_{x±} (*ϕ*
_{y±}) indicate the reflectivity and the total phase at *x*=±*L*/2 (*y*=±*L*/2), respectively. Applying the boundary conditions (11) to Eq. (6), we obtain the numerical solutions with the facet reflectivity. Note that, as in the previous section, we set the two coupling constants as κ_{2}L=-4, κ_{3}L=-2, and the surface-coupling coefficient as κ_{0}L=0.09.

Since these devices commonly have homogeneous material around the PC cavity, we assume the reflection *r* is uniform on all the facets. Figure 9(a) shows the threshold gain of mode A as a function of the reflectivity for the two phase conditions; *ϕ* =0 and *ϕ* =*π*. In the case of *ϕ* =π, the αL of mode A monotonically decreases with increasing ρ, while in the case of *ϕ*=0, the αL has a peak at around ρ=0.3. Figure 9(b) shows the threshold gain of the mode A as a function of the phase for the reflectivity of ρ=0.2. The curve peaks at around *ϕ*=1.9π. and takes it’s minimum at around *ϕ*=0.8π. This result provides a design rule for the facet phase to minimize the threshold gain of the 2D PC lasers.

## 5. Conclusion

We have developed the coupled wave theory for the 2D square lattice structure with TM polarization. Numerical calculations of the eigenvalue problem have shown the threshold gain, the wavenumber deviation and the mode pattern of the 2D resonant modes. We have classified the resonant modes by their longitudinal and transverse mode patterns and shown that the basic modes with the field distribution of the lowest order correspond to the resonant modes at the photonic band edge. The effect of the boundary reflector has also been presented for uniform phase and reflectivity. The results obtained in this paper provide fundamental insight into the 2D DFB effect of the PC lasers. A further development in designing and optimizing the 2D PC lasers by the current method is envisaged.

## Acknowledgments

This work was partly supported by Core Research for Evolutional Science and Technology – Japan Science and Technology Agency (CREST-JST). K. Sakai was supported by a Research Fellowship of the Japan Society for Promotion of Science.

## References and links

**1. **M. Imada, S. Noda, A. Chutinan, T. Tokuda, M. Murata, and G. Sasaki, “Coherent two-dimensional lasing action in surface-emitting laser with triangular-lattice photonic crystal structure,” Appl. Phys. Lett. **75**,316–318 (1999). [CrossRef]

**2. **M. Meier, A. Mekis, A. Dodabalapur, A. Timko, R. E. Slusher, J. D. Joannopoulos, and O. Nalamasu, “Laser action from two-dimensional distributed feedback in photonic crystals,” Appl. Phys. Lett. **74**,7–9 (1999). [CrossRef]

**3. **S. Noda, M. Yokoyama, M. Imada, A. Chutinan, and M. Mochizuki, “Polarization Mode Control of Two-Dimensional Photonic Crystal Laser by Unit Cell Structure Design,” Science **293**,1123–1125 (2001). [CrossRef] [PubMed]

**4. **E. Miyai, K. Sakai, T. Okano, W. Kunishi, D. Ohnishi, and S. Noda, “Lasers producing tailored beams,” Nature **441**,946 (2006). [CrossRef] [PubMed]

**5. **G.A. Turnbull, P. Andrew, W.L. Barnes, and I.D.W. Samuel, “Operating characteristics of a emiconducting polymer laser pumped by a microchip laser,” Appl. Phys. Lett. **82**,313”315 (2003). [CrossRef]

**6. **I. Vurgaftman and J. Meyer, “Design Optimization for High-Brightness Surface-Emitting Photonic-Crystal Distributed-Feedback Lasers,” IEEE J. Quantum. Electron. **39**,689–700 (2003). [CrossRef]

**7. **M. Imada, A. Chutinan, S. Noda, and M. Mochizuki, “Multidirectionally distributed feedback photonic rystal lasers”, Phys. Rev. B **65**, 1953061–8 (2002). [CrossRef]

**8. **K. Sakai, E. Miyai, T. Sakaguchi, D. Ohnishi, T. Okano, and S. Noda, “Lasing Band-Edge Identification for a Surface-Emitting Photonic Crystal Laser,” IEEE J. Sel. Areas Commun. **23**,1335–1340 (2005). [CrossRef]

**9. **M. Yokoyama and S. Noda, “Finite-Difference Time-Domain Simulation of Two-Dimensional Photonic Crystal Surface-Emitting Laser,” Opt. Express **13**,2869–2880 (2005). [CrossRef] [PubMed]

**10. **H. Kogelnik and C. V. Shank, “Coupled-Wave Theory of Distributed Feedback Lasers,” J. Appl. Phys. **43**,2327–2335 (1972). [CrossRef]

**11. **G. Scamarcio, F. Capasso, C. Sirtori, J. Faist, A. L. Hutchinson, D. L. Sivco, and A. Y. Cho, “High-Power infrared (8-Micrometer Wavelength) Superlattice Lasers,” Science **276**,773–776 (1997). [CrossRef] [PubMed]

**12. **R. Dorn, S. Quabis, and G. Leuchs, “Sharper focus for a radially polarized light beam,” Phys. Rev. Lett. **91**, 2339011–4 (2003). [CrossRef]

**13. **M. Plihal and A. A. Maradudin, “Photonic band structure of two-dimensional systems: The triangular lattice,” Phys. Rev. B **44**,8565–8571 (1991). [CrossRef]

**14. **H. Kogelnik, “Coupled Wave Theory for Thick Hologram Gratings,” Bell System Tech. J. **48**,2909–2947 (1969).

**15. **K. Sakai, E. Miyai, and S. Noda, “Coupled-wave model for square-lattice two-dimensional photonic crystal with transverse-electric-like mode,” Appl. Phys. Lett. **89**, 0211011–3 (2006). [CrossRef]

**16. **R. F. Kazarinov and C. H. Henry, “Second-order distributed-feedback lasers with mode selection provided by first-order radiation losses,” IEEE J. Quantum Electronics **21**,144–150 (1985). [CrossRef]

**17. **S. R. Chinn “Effects of Mirror Reflectivity in a Distributed-Feedback Laser,” IEEE J. Quantum Electron. **9**,574–580 (1973). [CrossRef]

**18. **W. Streifer, R. D. Burnham, and D. R. Scifres, “Effect of External Reflectors on Longitudinal Modes of Distributed Feedback Lasers,” IEEE J. Quantum Electron. **11**,154–161 (1975). [CrossRef]