1. Introduction

Photonic crystals (PCs) are periodic dielectric structures that have attracted much attention for enabling the control of radiation phenomena in wavelength-scale [1, 2]. Particularly, anomalously narrow resonances can occur in the continuum of PC slabs when the radiation is completely canceled by the destructive interference, which can be utilized in many applications such as filters [3, 4], lasers [5], and sensors [6]. Such a phenomenon has also been interpreted as a photonic analogy of the bound-states in the continuum (BICs) in quantum mechanics that was first proposed by von Neumann and Wigner [7] and was studied in decades [8, 9]. Recently, researchers have expended significant effort to construct BICs in photonic systems [10–15], and both at-Γ [16,17] and off-Γ BICs have been experimentally observed [18–20].

Owing to the protection of symmetries [17–22] or the cancellation of radiation at PC boundaries [23], the theoretical lifetime of a BIC is infinite in an ideal periodic PC slab. Unique high-Q resonances have been reported in large area or even macroscopic PC slabs [17]. However, since a realistic PC slab device always has a finite footprint and the structural fluctuations are inevitable due to the imperfection of the fabrication process, the Q factors of BICs degrade and become finite. As a result, it could be important to investigate the impact of structural fluctuations on the radiation properties of BICs, as well as its relationship with the footprint size of the PC slab. Conventionally, the statistical properties of Q factors were obtained by performing numerical simulations on a series of random distributed fluctuation patterns [24, 25], in which no doubt substantial computational resources are required.

In recent years, we developed a coupled-wave theory (CWT) framework to realize accurate and efficient analysis of the detailed wave interactions within the PC slab [26–28]. The framework is capable of analyzing finite [29] or infinite PC slabs with different crystalline geometries and arbitrary tilted sidewalls [27] on both TE- and TM-like modes [30–33]. Moreover, analytical expressions of radiative waves can be obtained to derive the condition of BICs [23]. In this work, we present an analytical and statistical investigation on structural fluctuations induced radiation in PC slabs, based on the CWT framework.

From the perspective of CWT, the radiations caused by structural fluctuations can be treated as a perturbation to the BICs. For a perfect periodic structure, the permittivity as well as the electromagnetic field can be expanded into a series of Bloch functions in the basis of integer multiple of the primitive reciprocal lattice vector. The existence of structural fluctuations would break the periodic condition of a unit cell, contributing to some fractional orders of Bloch waves. Some of the fractional order waves are radiative and will transport energy away from the slab. As a result, the Q factors of BICs degrade.

To investigate the radiation induced by the structural fluctuations, we first extend the CWT framework to a supercell lattice which contains N periods of unit cells. The permittivity within a supercell lattice is randomly perturbed from the ideal periodicity according to normal distributions with a given standard deviation. The structural fluctuations are depicted by such a standard deviation which can be determined from the realistic fabrication precision. Furthermore, the fractional order radiative waves can be calculated for analytical and statistical evaluations of the Q factors under perturbations. When N becomes large enough, the radiation of a PC slab becomes dominated by structural fluctuations, whereas the impact of the periodic condition turns to be negligible, and consequently the supercell lattice model approaches a realistic finite size PC device with structural fluctuations.

The remainder of this paper is organized as follows. In Section 2, we derive the CWT formulation within a supercell lattice. In Section 3 and 4, we discuss the coupling paths and the radiation loss by assuming independently normally distributed structural fluctuations. In Section 5, we present our results and discussions. In Section 6, we conclude with our findings.

2. Theory and formulation

A typical one dimensional (1D) PC slab is shown in Fig. 1. It is a three layered structure, periodic in the x direction, and extends infinitely in the y direction. The structural parameters we use in this work is presented in Table 1. Since an ideal structure is strictly periodic in the x direction, the electromagnetic field can be expanded as the sum of a series of Bloch waves. For TE modes, the electric field is only polarized in the y direction, namely E = Ee_y, and E follows $E(x,z)={\sum }_m{E}_m(z)e^{-im{\beta }_{0}x}$, m ∈ Z.

 

Fig. 1 (a)(b) The schematics of a typical 1D-PC slab. It is periodic in the x direction, extends infinitely in the y direction, and is three layered in the z direction. The permittivity of the upper/lower cladding layers are ε_u and ε_l, respectively. The PC layer is patterned with two materials ε_a and ε_b with a period of a and a filling factor of f₀ = w/a. (c) Two types of structural fluctuations concerned in this work: (i) fluctuations of the center positions of the holes, (ii) fluctuations of the filling factors. (d) A typical band structure near the 2nd order Γ point matched with mode TE₀. We denote the modes as A and B in increasing frequency.

Download Full Size | PDF

Table 1. Structural Parameters

View Table

For a realistic structure, structural fluctuations are inevitable due to the imperfection of the fabrication process, which breaks the ideal periodicity. Here we consider two types of structural fluctuations: (i) fluctuations of the center positions of the holes, (ii) fluctuations of the filling factors. The fluctuations follow normal distributions with given standard deviations, which can be determined by the realistic fabrication precision.

Since the ideal periodicity is broken, the Bloch theorem doesn’t apply any more. Here, we propose a supercell lattice extension of the CWT framework, to deal with the structural fluctuations. We assume that the structure is still periodic in a larger scale, which we refer to as a supercell. Structural fluctuations happen in arbitrary ways within a supercell, and the supercell repeats itself periodically in the whole infinite structure. Therefore, the Bloch Theorem still holds for the supercell but fractional orders of Bloch waves occur. In this case, the electric field is expressed as $E(x,z)={\sum }_m{E}_m(z)e^{-im{\beta }_{0}x}$, m ∈ Z/N, where N is the size of a supercell in the number of primitive lattice constant. Meanwhile, the m-th order Fourier coefficient of the permittivity ε in the PC layer is given by ε_m = ε_z,m + Δε_m, m ∈ Z/N, where ε_z,m and Δε are the Fourier transformation of the ideal ε without fluctuations, and of the structural fluctuations itself, respectively, given by:

The modes in the PC slab are phase-matched waveguide modes. The phase-matching condition can be fulfilled on any order of the waveguide mode (TE₀, TE₁, … and TM₀, TM₁, …) and on any order of the Γ point. In this work, we focus on the 2nd order Γ point and the TE₀ mode.

Therefore, in the unitcell CWT framework, the Bloch waves of orders m = ±1 occupy most of the whole energy near the 2nd order Γ point, and they are referred to as basic waves. We use them as the source to excite Bloch waves of other orders, and at last obtain the eigen frequencies and eigen vectors by solving a eigenvalue problem of a coupling matrix. The eigen value presents a complex mode frequency, and hence, we obtain the Q factor Q = f_r/2/f_i, where f_r and f_i are the real and imaginary parts of the eigen frequency, respectively.

In the extended supercell CWT, as N grows large, the k-space distance of adjacent fractional Bloch waves 1/N approaches 0, and the fractional order waves near m = ±1 become quite close to the phase-matching condition of β = β₀. As a result, the coupling strengths between them become large, and they also occupy sufficient large energy. In this case, we need to adopt the degenerate perturbation method rather than previous used non-degenerate perturbation method. More specifically, Bloch waves of m ≈ ±1 instead of only that of m = ±1 should be chosen as the basic waves. Still, we assume that the basic waves have the same profile as the phase-matched waveguide mode, i.e., E_m(z) = V_mΘ_m(z), where m ∈ V, V is the set of the basic waves, and Θ_m is the profile of the phase-matched waveguide mode of order m. Θ_m(z) satisfies:

We calculate Eq. (5) − Eq. (4)×V_m, and multiply it by ${\mathrm{\Theta }}_m^{*}(z)$, and then integrate it along z ∈ (−∞, ∞), then we obtain:

Equation (7) leads to an eigenvalue problem. Since the wavenumber k on the right hand side of Eq. (7) affects the coupling strength, it needs to be evaluated in advance and we choose it as the wavenumber without fluctuations, referred to as k_c. For mode A (high-Q) we are concerned with, k_c is the wavenumber of mode A without fluctuations, and can be obtained from the unicell CWT.

Let matrix $D_{m{m}^{\prime }}=k_m^2{\delta }_{m{m}^{\prime }}$ and $L_{V,m{m}^{\prime }}=-(k_c^2/h_m){\epsilon }_{m-m^{\prime }}{O}_{m{m}^{\prime }}(1-{\delta }_{m{m}^{\prime }})$, where m, m′ ∈ V; $L_{A,m{m}^{\prime }}=-(k_c^2/h_m){\epsilon }_{m-m^{\prime }}{O}_{m{m}^{\prime }}$, where m V, m′ ∈ A, Eq. (7) is rewritten in a matrix form:

Further, we can use Green’s function method to obtain the column vector A as a multiplication of a matrix with V. Define the Green’s function

Solving this eigenvalue problem, we obtain the complex wavenumbers and the corresponding eigen vectors, and hence the mode frequencies and Q factors. Moreover, for the Bloch waves of m ∈ A, their amplitudes can be obtained by Eq. (13), with their profiles following:

From the perspective of quantum mechanics, the energy levels with close energies should better be dealt with the degenerate perturbation method instead of the non-degenerate perturbation method. In the context of CWT, it means the Bloch waves of m ≈ ±1 should be included in the V. The selection of the basic waves (set V) depends on the coupling strength from the Bloch waves of m = ±1 to their neighbors. In the computations we implement in this work, we choose Bloch waves of orders 0.95 ≤ |m| ≤ 1.05 into set V, unless specifically declared.

It is noteworthy that comparing with our previous CWT, one significant improvement in this work is that we can analytically deal with the iterations. As elaborated previously [28, 32], owing to the perturbation method nature of CWT, it requires some iteration techniques to obtain quantitative accurate self-consistent solutions when the structure has high-index-contrast. In this work, we substituting the field iteration with the inversion of a matrix. This provides two advantages:

First, it overcomes the issue of bad convergence of iterations in some situations. More specifically, the iteration is equivalent to expanding (1 − C)⁻¹ in Eq. (13) to 1 + C + C² + C³ + …. Good and efficient convergence requires that the modulus of every eigenvalue of matrix C is far less than 1, which cannot be always fulfilled in the supercell case. Iteration doesn’t converge for some random fluctuation patterns. In contrast, the matrix inversion operation in Eq. (13) always holds as long as the determinant of (I − C) is not zero.

Second, matrix inversion operation takes only one matrix operation and can be accelerated by mature computational libraries. In contrast, the iteration procedure needs to renew all the field profiles and re-calculate the matrix every time till convergence. As a result, the inversion of the matrix is far more fast and efficient than performing iterations. As a comparison, for a supercell structure of N = 100, the proposed method costs 5.49s, whereas the iteration procedure costs about half an hour for 20 iterations to obtain the result.

3. Coupling paths in a fluctuated structure

When applying Eq. (14) to computation, a proper truncation of set A has to be made. That is, we only take the Bloch waves of orders |m| ≤ M. Our computation shows that M = 5 is enough for the structure we work on. With such a truncation, the computational complexity is O(N²M²). However, with the structural fluctuations presented in Eq. (2), some approximations can be adopted to optimize the computational cost to O(n_V N M²), where n_V is the number of elements in set V. The main computational cost concentrates on the calculation of matrix T in Eq. (13). We first rewrite matrix T as:

Let C = C_z + C_d, C₀ = C_0z + C_0d, where

Equation (16) is rewritten as

Figure 2 illustrates the major coupling paths according to Eq. (21). For simplicity, It only shows the case when V = {1, −1}. Let U = A ∪ V represent all Bloch waves in the system. The elements in set U are classified into N groups according to their orders m; that is, group U⁽ⁱ⁾ contains the waves that satisfy (m mod N = i). Sets A and V are grouped according to the same rule, as A⁽ⁱ⁾ and V⁽ⁱ⁾.

 

Fig. 2 A schematic of coupling paths indicated in Eq. (21), by assuming the basic wave set V = {−1, 1}. The big dots represent the Bloch waves of m ∈ Z, whereas the small ones depict the fractional Bloch waves of m ∉ Z. The coupling paths indicated by the coupling matrices C_0z, C_0d, C_z, C_d are illustrated in the figure.

Download Full Size | PDF

The coupling paths can be divided into strong coupling paths (C_z and C_0z) which is dominated by ideal permittivity periodicity, and weak coupling paths (C_d and C_0d) which is caused by structural fluctuations. Strong coupling only happens inside a certain group U⁽ⁱ⁾. In contrast, weak coupling occurs between arbitrary two Bloch waves, namely, it induces the coupling between U⁽ⁱ⁾ and U^(j), where i ≠ j.

According to the analysis above, Eq. (21) can be calculated by groups; that is:

According to Fig. 2 and Eq. (22), the computation could be performed by each group first, and then gather them together. Denoting the number of elements in set V as n_V, the optimized computational complexity turns to be O(n_V N M²), which is much more efficient than O(N²M²). When applying the proposed grouping algorithm to two dimensional (2D) PC slabs, the complexity can be optimized from O(N⁴M⁴) to $O(n_V^2{N}^2{M}^4)$, which would remarkably reduce the computational resource.

4. Fluctuation induced radiation loss

The structural fluctuations will induce fractional orders of Bloch waves and some of them are radiative leaky waves. Figure 3 presents the amplitude distribution in the k-space of mode A on SOI under a random fluctuation pattern with supercell size N = 300. It is noticed that because the fluctuations have an average value of 0, the m = 0 order wave still doesn’t radiate, as in the ideal periodic structure. The radiation is mainly induced by fractional order Bloch waves, particularly near m = 0. Calculation shows that Bloch waves of −0.05 ≤ m ≤ 0.05 contribute over 90% of the whole radiation power. As mentioned above, the Bloch waves of m ≈ ±1 have considerable energy because they are almost degenerate with that of m = ±1. According to the analysis in Section 3, they excite the radiative waves of m ≈ 0 through strong coupling paths, whereas the radiative waves far away from m = 0 are mainly excited by basic waves though weak couping paths.

 

Fig. 3 The amplitude distribution in the k-space for mode A on SOI under a random fluctuation pattern, with supercell size N = 300, (a) near m = 0 (b) near m = 1, (c) near m = −1, and (d) large scale. The red area in (d) represents the Bloch waves above the light cone. Data is obtained from CWT.

Download Full Size | PDF

Since CWT gives an analytical expression of the radiative waves and the structural fluctuations are modeled as normal distributions with a given standard deviation, we can calculate the expected value of the radiation loss directly, given by 1/Q = 2 f_i/f_r. Because the fluctuated structure is treated as a perturbation of an ideal periodic structure, we can assume the real part f_r remains constant and the problem comes to calculate the expected value of f_i.

For a given set V, f_i can be evaluated from the leakage of energy as f_i = P/4πW, where W = (1/2)Σ_m(|V_m|² ∫ ε₀(z)|Θ_m(z)|² dz) is the modal energy in a unit cell, and P =Σ_m P_m is the radiation power. More specifically, P_m = P_m,u + P_m,l, where u, l represent the radiation into the upper and lower claddings, respectively. $P_{m,u}=\frac{1}{2}\sqrt{{\epsilon }_{0,u}}|E_{m,u}{|}^2\cos {\theta }_{m,u}$, and $P_{m,l}=\frac{1}{2}\sqrt{{\epsilon }_{0,l}}|E_{m,l}{|}^2\cos {\theta }_{m,l}$, where ε₀,_u and ε₀,_l are the permittivity of the upper and lower cladding respectively, θ_m,u and θ_m,l are the angles of the emergent radiation waves correspondingly. In addition, the radiation power of evanescent waves is zero. Noticing E_m,u = A_mΘ_m,u and E_m,l = A_mΘ_m,l, we have:

The overall expected radiation power is E(P) = Σ_m E(P_m). To obtain E(P) and E(f_i), we need to calculate E(|A_m|²). For mode A at the 2nd-Γ point, the first term in Eq. (21) doesn’t contribute to radiation, and hence Eq. (13) is rewritten as follows to calculate the radiation power:

Each element in C_0d and C_d is a perturbation item of the permittivity, and they can be expressed as: $C_{0d,ij}={\displaystyle \sum _l{a}_{ij,l}\mathrm{\Delta }{\rho }_l+}{b}_{ij,l}\mathrm{\Delta }{\xi }_l$, and $C_{d,ij}={\displaystyle \sum _l{c}_{ij,l}\mathrm{\Delta }{\rho }_l+}{d}_{ij,l}\mathrm{\Delta }{\xi }_l$, where the elements of matrices a, b, c, and d are derived in Appendix. Denoting V_z = C_0z V, Eq. (24) has the following form as derived in Appendix:

Here, we choose the Bloch waves of m = ±1 as basic waves, i.e., V = {1, −1}. In addition, we assume that Δξ_l has an average value of 0; that is, define $\mathrm{\Delta }{\xi }_l=\mathrm{\Delta }{\xi }_l^{(0)}-\mathrm{\Delta }{\overline{\xi }}^{(0)}$, where $\mathrm{\Delta }{\xi }_l^{(0)}$ = Δx_l/a, Δx_l is the variation of the center positions of the holes under structural fluctuations. This definition translates the PC layer by $(-{\mathrm{\Delta }}_{\xi }^{-(0)})$ along the x-axis. It affects nothing except the phase of V. However, under this definition, the eigen vector V remains almost invariant during perturbation, and the statistical terms only lie in matrices C_0d and C_d. This makes the calculation of the expected value much easier.

By assuming independently normal distribution of Δρ_l and $\mathrm{\Delta }{\xi }_l^{(0)}$, namely $\mathrm{\Delta }{\rho }_l~N(0,{\sigma }_{\rho }^2)$, and $\mathrm{\Delta }{\xi }_l^{(0)}~N(0,{\sigma }_{\xi }^2)$, We obtain the following equation:

In this way, the expected value of 1/Q is obtained. Since h_il and g_il are independent on the fluctuation pattern, we find that the expected radiation power is proportional to ${\sigma }_{\rho }^2$ σ² and ${\sigma }_{\xi }^2$ from Eq. (26).

5. Results and discussions

According to the derivation above, we investigate the radiation properties of realistic 1D-PCs on SOI wafers with parameters listed in Table 1. We assume the fabrication process precision is about 1nm, which gives the standard deviation of the structural fluctuation σ = 1/550. For a given fluctuation pattern, the Q factor can be calculated from CWT analytically, or from numerical simulation of the finite-element method (FEM, COMSOL Multiphysics). On the other hand, the statistical properties of the Q factor under normally distributed random fluctuations can be evaluated from a Monte Carlo simulation by using both CWT and FEM. With a given type of fluctuation and a known standard deviation, CWT also gives a close form of the expected value of the Q factor as shown in Eqs. (23)–(26).

Figure 4 illustrates 1/Q of mode A under 25 randomly generated fluctuation patterns on an SOI supercell of N = 20, in the cases of varying Δρ and Δξ separately. For each pattern, the Q factors calculated by CWT and FEM agree well with each other, which strongly validates the proposed CWT method. Meanwhile, the expected value of 1/Q calculated by CWT shows good consistency with the average value of FEM simulation as well, indicating that the expected Q factor provides a quick and effective estimation of the radiation properties under structural fluctuations, without running tedious Monte Carlo simulation.

 

Fig. 4 1/Q of mode A on an SOI supercell of N = 20. As a Monte Carlo simulation, 25 randomly generated fluctuation patterns are plotted with (a) σ_ρ = 1/550, σ_ξ = 0, and (b) σ_ρ = 0, σ_ξ = 1/550. The results of both CWT (blue circles) and FEM (red plus sign) are presented. σ_ρ and σ_ξ are determined by assuming the fabrication process precision of 1nm. The green dashed lines show the expected value of 1/Q calculated by CWT, which agrees well with the average value of FEM results.

Download Full Size | PDF

For a certain pattern, we denote the fluctuations of filling factors and the center positions under normal distributions of a standard deviation σ₀ as $\mathrm{\Delta }{\rho }_l^{(0)}$ and $\mathrm{\Delta }{\xi }_l^{(0)}$, respectively. For a different standard deviation σ, the fluctuations turn to be Δρ_l and Δξ_l; we have $\mathrm{\Delta }{\rho }_l=\mathrm{\Delta }{\rho }_l^{(0)}\cdot \sigma /{\sigma }_0$ and $\mathrm{\Delta }{\xi }_l=\mathrm{\Delta }{\xi }_l^{(0)}\cdot \sigma /{\sigma }_0$. The Q factors versus increasing σ_ρ and σ_ξ are presented in Fig. 5, in which a quadratic dependence of 1/Q on a varying σ is illustrated. Such results agree well with the prediction of CWT. For a given pattern, we derive from Eq. (25) as:

 

Fig. 5 Dependence of 1/Q of mode A on σ² under five random fluctuation patterns on an SOI supercell structure of N = 20. Results of both CWT and FEM are presented when fluctuations exist in (a) only Δρ, (b) only Δξ, and (c) both Δρ and Δξ. Each of them shows a clear quadratic dependence of 1/Q on σ.

Download Full Size | PDF

The above calculations have assumed a relatively small supercell size (N = 20). However, a realistic structure always contains more periods. Figure 6 shows the radiation loss of mode A on SOI supercell structure, under random fluctuated patterns of increasing N. Again, the results of CWT and FEM agree well with each other. It is noteworthy that, a degenerate perturbation of basic wave set V = {m | 0.95 ≤ |m| ≤ 1.05} CWT is adopted for large N. In this case, the expected radiation loss Eq. (26) depends on summing up a large number of random variables, and thus becomes less obvious. Therefore, we alternatively calculate the expected value of 1/Q by running and averaging a series of Monte Carlo simulations of CWT. As the shown in Fig. 6, the expected radiation loss 1/Q approaches a stable asymptotic value when N goes large, indicating that further increase the supercell size will not promote the Q factor any longer.

 

Fig. 6 1/Q of mode A on an SOI supercell structure under random fluctuation patterns versus varying supercell size N, with σ_ρ = 1/550 and σ_ξ = 1/550. Results calculated by both CWT (blue circles) and FEM (red plus signs) are presented for individual patterns. The gray dashed line indicates the critical supercell size N = 67 of switching non-degenerate perturbation to degenerate perturbation. The purple line shows the asymptotic value of 1/Q.

Download Full Size | PDF

The criterion of choosing degenerate perturbation or non-degenerate perturbation relies on whether or not the fractional order Bloch waves m = ±1 occupy a large portion of the modal energy. More specifically, among the Bloch waves 0.9 ≤ |m| ≤ 1.1, if the Bloch waves of m = ±1 occupy less than 90% of the total energy, we turn to adopt the degenerate perturbation. Under fluctuations of σ_ρ = σ_ξ = 1/550 on an SOI wafer, it gives a critical supercell size of N = 67.

When N grows large enough, we expect the radiation loss turns to be dominated by the structural fluctuations, and the influence of the boundary condition becomes considerably weak. In this case, there turns to be no difference in adopting the periodic boundary condition or the finite-size boundary condition. As a result, the supercell lattice model approaches a realistic finite size PC device. Figure 7(a) presents the radiation loss of mode A under 25 random fluctuation patterns by adopting supercell structure of N = 300, and the corresponding finite-size structure with the same fluctuation patterns. The finite-size structure and the supercell structure with the same size and fluctuation patterns give the same radiation losses, which supports that the proposed supercell CWT model is capable of dealing with finite-size structures with large N.

 

Fig. 7 (a) The CWT and FEM calculated 1/Q of mode A on SOI supercell structure under 25 random fluctuation patterns, with supercell size N = 300, and the FEM results for the corresponding finite-size structure with the same fluctuation patterns. (b) The asymptotic 1/Q versus the structural fluctuation σ obtained from CWT using Monte-Carlo method. Data is obtained by setting the supercell size N = 300.

Download Full Size | PDF

Figure 7(b) shows the asymptotic value of radiation loss under different fabrication precisions. According to the analysis above, the Q factors of a large-area PC device is limited by the fabrication precision, i.e. the standard deviation of the structural fluctuations. When the device size N exceeds certain periods, increasing N would no longer benefit the Q factors since the fluctuation-induced scattering loss dominates. Such phenomenon was recently reported in periodic array of dielectric spheres [15]. For the SOI PC structure we focus on, 300 periods (Fig. 6) are large enough to ensure the Q factor to approach its asymptotic value under fabrication precisions from 0.5nm to 1nm. For typical fabrication precision of 1nm, the Q factor degrades from almost infinite to ∼ 3 × 10⁴.

6. Conclusion

In this work, we extend the CWT framework to a supercell lattice PC structure in order to model structural fluctuations which are inevitable in realistic devices. The permittivity distributions are randomly perturbed from ideal periodicity due to the imperfection of fabrication process, and hence leading to the existence of fractional orders of Bloch waves. The extended CWT model is capable of depicting the details of wave interactions within the fluctuated PC structure, allowing us to perform analytical and statistical evaluations on the radiation properties of high-Q resonances.

Because some of the fractional order waves lie above the light-cone, they are radiative and will transport energy away from the slab. The high-Q resonances which are supposed to have infinite lifetime in ideal periodic structure degrade to lower Q factors. For a given fabrication precision, we investigate two type of structural fluctuations, i.e. (i) fluctuations of the center positions of the holes, (ii) fluctuations of the filling factors, by assuming that the structural parameters obey normal distributions with given standard deviations. Accordingly, the coupling paths as well as the expected Q factors are discussed.

For given fluctuated patterns, the Q factors calculated from CWT show good consistency with the FEM simulations. Moreover, the statistical properties of the radiation loss are evaluated from running a series of Monte Carlo simulation by using CWT and FEM respectively. The average values of Monte Carlo simulation agree well with the expected values from CWT. These results confirm the validity of the supercell CWT. Besides, we notice that, the difference between periodic boundary condition and finite-size boundary condition turns to be negligible for a large enough supercell size, and consequently the supercell lattice model approaches a realistic finite area PC with structural fluctuations.

From the perspective of CWT, we find that the Q factors of high-Q resonances are constraint by the structural fluctuations in realistic PC devices, and they will converge to a certain value when the devices have large enough periods. For a typical fabrication precision of 1nm on a 300 periods SOI PC slab, the Q factor degrades to 3 × 10⁴. Further increasing the device size will no longer promote the Q factor. As a result, to realize high-Q resonances in PC structures for filter, laser and sensor applications, suppressing the structural fluctuations during fabrication is essentially important.