## Abstract

The perturbation method is combined with the Rigorous Coupled Wave Analysis (RCWA) to enhance its computational speed. In the original RCWA, a grating is approximated by a stack of lamellar gratings and the number of eigenvalue systems to be solved is equal to the number of subgratings. The perturbation method allows to derive the eigensolutions in many layers from the computed eigensolutions of a reference layer provided that the optical and geometrical parameters of these layers differ only slightly. A trapezoidal grating is considered to evaluate the performance of the method.

© 2010 Optical Society of America

## 1. Introduction

For many years the electromagnetic wave diffraction by surface relief gratings has been widely studied. Most of the methods which deal with this subject are based on a differential approach and use the staircase approximation of the grating profile. This technique of considering the grating as a stack of lamellar gratings was first suggested by Peng *et al.* [1] and was later used in the very popular Rigorous Coupled Wave Analysis (RCWA) [2]. In all the methods employing this technique, the electromagnetic field is written as a modal expansion in each layer and the expansions in the different layers are connected by boundary conditions. Numerical difficulties which may be encountered are linked, on the one hand, to the way in which the field is calculated inside the lamellar gratings, and, on the other hand, to the way in which the field expansions are propagated through the layers.

The field components in a lamellar grating may be thought of in terms of modal field expansions. The classical modal method [3, 4] is a rigorous theory since it uses basis functions which constitute exact solutions of Maxwell’s equations and satisfy the pseudo-periodic boundary conditions. The Fourier modal method (FMM) was introduced by Knop [5] in the case of dielectric gratings. In this technique, each mode is approximated by a truncated Fourier series. The application of this method to a metallic grating in TM polarization was much questioned until 1996 when several authors [6–8] proposed a dramatic improvement based on an efficient Fourier expansion of the periodic permittivity function. In contrast with the classical modal method, the FMM leads to a numerical scheme which permits straightforward implementation.

This feature explains the popularity of the RCWA which combines the multilayer approach with the FMM, whereas the numerical implementation of the classical modal method requires more sophisticated mathematics in order to overcome numerical difficulties [9,10]. In all cases, convergence with respect to the number of layers depends on the way in which the modal field expansion is propagated through the layers [8, 11, 12]. Many studies have been devoted to this subject up to now [**?**]. In contrast to the so-called differential method [14] and the coordinate transformation method (C-method ) [15–17], all the above mentioned methods rely on the assumption that the staircase approximation is valid. It has been shown that it is not the case for highly conducting gratings. Although the RCWA fails in such cases, it remains very efficient in most cases thanks to both its versatility and simplicity.

The staircase approximation entails solving as many eigenvalue problems as there are layers. The aim of this paper is to show that the number of eigenvalue problems can be drastically reduced by means of a perturbation method (PM). This approach is based on the fact that the permittivity function varies very little from one layer to the adjacent one.

## 2. Staircase approximation

Let us consider the grating depicted in Fig. 1. This structure is invariant along the *z*-direction and it is illuminated by an incident wave, under incidence angle *θ*, that can be either TE (the only non nul components of the field are *E _{z}*,

*H*and

_{x}*H*) or TM (the only non nul components of the field are

_{y}*H*,

_{z}*E*and

_{x}*E*) polarized. It is well known that in these cases of polarization, all the field components can be expressed in terms of the

_{y}*E*-component in the TE case or in terms of the

_{z}*H*-component in the TM case. Let us add that throughout this paper, the time dependence of the fields will be held by the term exp(

_{z}*iωt*) where

*ω*denotes the circular frequency of the monochromatic incident wave. The incident medium is vacuum and the substrate is a dielectric medium with or without losses.

*E*and

_{z}*H*satisfy the propagation equations:

_{z}*k*= 2

*π/λ*denotes the vacuum wave number and

*ν*is the refractive index.

The grating surface is approximated by a staircase profile (see Fig. 2) so that the structure is replaced by a stack of layers, each of them is characterized by its thickness *h _{j}* and its refractive index

*ν*which only depends on the

_{j}*x*variable:

## 3. Eigenvalue problem

In this section, we solve the propagation equation Eq. (1) in any *j ^{th}* layer, therefore the subscript index

*j*is omitted. The equations Eq. (1) can be rewritten in a concise form:

*F*(

*x,y*) denotes

*E*or

_{z}*H*. ℒ is an operator which is different according to the polarization:

_{z}Since the operator ℒ depends only on the *x*-variable, the separation of variables method leads to

*ψ*is the eigenfunction of the ℒ-operator: that satisfies the pseudo-periodic condition:

*α*is determined by the incident wave:

*α*= sin

*θ*. The square root of the eigenvalue

*r*

^{2}is chosen in such a way that

Let us recall briefly some important properties of the ℒ-operator [3, 4] that will be useful for the perturbation scheme that will be developed in section (5). For that purpose, we define an inner product for any pseudo-periodic functions *f* and *g*:

*f*denotes the complex conjugate of

*f*.

By using the pseudo-periodicity, it can be established that the operator ℒ^{†} defined by:

*ψ*and ${\psi}_{m}^{\u2020}$ are orthogonal with respect to the convenient inner product. The completeness of eigenfunctions allows to expand any pseudoperiodic function

_{m}*f*(

*x*) as follows:

*F*(

*x,y*) (

*i.e. E*or

_{z}*H*) may be expanded in the eigenfunction series:

_{z}According to Eq. (8)*A _{m}* is associated with downgoing waves and

*B*to upgoing waves. When

_{m}*ν*(

*x*) is a real function ℒ

^{†}and ℒ are identical i.e. ℒ is self-adjoint.

## 4. Fourier modal expansion

The eigenfunctions can be calculated exactly by solving a transcendental equation [3, 4]. Another way for implementing the modal method [2, 5] consists in expanding any eigenfunction *ψ _{m}* in the Fourier series:

Numerically only a finite sequence of *N* Fourier coefficients is kept. Equation (6) involves products of functions of the *x*-variables which are represented by truncated Fourier series. Li described the procedures for correctly converting such products in discrete Fourier space [8]. Solving Eq. (6) leads to the research of the eigenvalues and eigenfunctions of a (2*N* + 1) × (2*N* + 1) matrix:

The *ε _{mn}* and

*ε̃*elements of ℰ and ℰ̃ respectively are the (

_{mn}*m*–

*n*)

*Fourier coefficients of*

^{th}*ν*

^{2}(

*x*) and 1

*/ν*

^{2}(

*x*). Each eigenvalue of

**L**

*and*

_{E}**L**

*is associated with an approximate*

_{H}*r*eigenvalue which is computed according to Eq. (8). Finally, the modal expansion of Eq. (5) is approximated as follows:

_{m}*ψ*denotes

_{nm}*n*Fourier coefficients of the eigenfunction

^{th}*ψ*(

_{m}*x*).

## 5. Perturbations

Perturbation theory is a mathematical method used to find an approximate solution to a problem from the exact solution of a related problem. This method is applicable if the problem at hand can be formulated by adding a “small” term to the mathematical description of the exactly solvable problem. Historically, the perturbation method has its roots in early celestial mechanics and it was first used to solve algebraic equations, before being applied to the operator theory, especially in classical quantum mechanics [19, 20].

Assuming that the eigenvalues and eigenfunctions of Eq. (6) are exactly known, we introduce the superscript index (0) to refer to this case which is characterized by the refractive index function *ν*^{(0)}(*x*). In accordance with this notation Eq. (6) becomes:

*p*. In any layer, the refractive index function

*ν*(

*x*) is slightly different from

*ν*

^{(0)}(

*x*), thus the operator ℒ may be expressed using the unperturbated operator ℒ

_{0}plus a correction term:

*ζ*is presumed to be “small”. It must be emphasized that the only operator

*P*is known and

*ζ*is a parameter useful only for analytical calculus. The eigenvalues ${\lambda}_{p}={r}_{p}^{2}$ and eigen-functions

*ψ*of the operator are expressed in terms of power series of

_{p}*ζ*:

*order approximation ${\lambda}_{p}={\lambda}_{p}^{(0)}$, ${\psi}_{p}={\psi}_{p}^{(0)}$.*

^{th}To calculate the higher order perturbation correction terms, we put the equations Eq. (19) and Eq. (20) into Eq. (6):

By successively identifying the 0, 1 and 2 rank coefficients of *ζ* we get:

Afterwards we consider the case where the operator is self-adjoint which corresponds to a real *ν*(*x*) index. The eigenvalues are therefore real. The complex eigenvectors are normalized by setting conditions on both the modulus and the phase:

By applying this normalization to 0, 1 and 2 orders we get:

In order to arrive at the first order approximation the equation Eq. (22b) is projected on ${\psi}_{p}^{(0)}$. By using Eq. (24a) and Eq. (24b) as well as the self adjoint property of the operator we get the expressions of the eigenvalues

and the eigenvectorsBy combining Eq. (22c) and Eq. (24c) we obtain the second order approximation of the eigenvalues

## 6. Numerical results

In this section, the perturbation method (PM) is applied to a surface grating made of a stack of two dielectric trapezoidal elements. This grating is characterized by six geometrical parameters indicated on Fig. (3): *d* = 480 nm, *b*_{1} = 240 nm, *b*_{2} = 120 nm, *b*_{3} = 80 nm, *h*_{1} = 80 nm, *h*_{2} = 200 nm. The direction of the incident plane wave is given by *θ* = 75° and the refractive index of the grating is *ν* = *n*_{1} = 1.51309 for the wavelength *λ* = 240 nm.

First, we investigate the convergence of the RCWA with respect to the truncation order *M* (Table 1) and the number of layers *N _{c}* (Table 2) for the zero-order reflected efficiency. These results show that, for the studied case, the convergence of the RCWA is slow, with respect to

*M*and to

*N*, in TE polarization as well as in TM polarization. In order to assure an accuracy equal to 10

_{c}^{−3}for the considered case, we will use throughout this paper the following parameters:

*Nc*= 220 and

*M*= 30.

Secondly, these results are taken as reference for evaluating the performance of the perturbation method. This is why numerical results are given with 15 digits even if only 3 digits are significant. In the perturbation method, the number of layers *Nc* is preserved but the eigenvalue problem is only solved for *n _{p}*-number of layers smaller than

*N*. The solutions for the remaining layers are obtained from Eqs. (25)–(28) at different perturbation orders. We define an error function

_{c}*ξ*as:

Here Int(*x*) provides the integer part of *x*, eff_{RCWA} and eff_{PM} refer to the zero-order reflected efficiency with RCWA and perturbation method respectively. We must note that the solving procedure of equations deduced from boundary conditions is the same as the one intervening in the RCWA. The first column of Tables 3–8 list the number of equations *n _{p}* that have been solved and the last column presents the time-saving

*ts*(

*n*) = 1 –

_{p}*t*

_{PM}(

*n*)

_{p}*/t*

_{RCWA}.

*t*

_{PM}is the computation time of the PM and

*t*

_{RCWA}is referred to the RCWA method.

The results obtained by approximating the eigenvalues and eigenvectors at first order, by means of Eq. (25) and Eq. (26), are presented in Tables 3 and 4. We notice that the 10^{−3} accuracy,as defined by Eq. (29) is obtained for a value of *n _{p}* = 4 i.e. by solving only 4 eigenvalue problems instead of 220 for both polarizations. Let us emphasize that the 10

^{−3}accuracy is the one of RCWA reference method. The corresponding time-saving amounts to 25%. In addition, we can study the intrinsic efficiency of the perturbation method in relation to the RCWA method by supposing that the results provided by the latter method are accurate. We can remark that the PM converges with respect to

*n*and that 10

_{p}^{−4}accuracy is achieved for

*n*= 20 and with a time-saving amounting to 25%.

_{p}The results obtained with the second-order approximation of both he eigenvalues ans eigenvectors i.e. by means of Eq. (27) and Eq. (28) are given in Tables 5 and 6. These results show that the second-order PM does not provide any significant improvement over the first-order one. For the TE polarization, we get 10^{−3} accuracy by solving fewer eigenvalue equations than the first order (2 instead of 4) but with a time-saving of only 13% instead of 25%. In others words, second order approximation is formally better than first-order one but this amelioration is numerically annihilated by the complexity of the second order approximation formulas. The situation worsens for TM polarization. By comparing the results in Table 6 with those in Table 4, we can notice that the second order approximation does not provide a better accuracy than first order one. The comparison of the operators in TE polarization [Eq. (15a)] and TM [Eq. (15b)] can explain these results. Only the *l _{j}* parameter (Cf. Fig. 2) differs from one layer to another. This parameter has influence through only the

*ν*

^{2}(

*x*)-function on the TE operator while it has influence more through the function 1/

*ν*

^{2}(

*x*) on the TM operator. Equations (26) and (28) show that the computation of the corrective term of second order eigenvectors increases the computation time. This observation has led us to numerically test a method that uses the correction of eigenvalue at second order Eq. (27) combined with the correction of eigenvectors with first order Eq. (26). These results are shown in Tables 7 et 8. Overall the results do not differ much from those obtained with the first order method as far as accuracy and computation time are concerned.

## 7. Conclusion

We have shown that the perturbation method allows to improve the run time of codes based on RCWA. Indeed, in many layers modal expansions may be deduced from those calculated in a reference layer thus avoiding to solve eigenvalue problems which is the most consuming part of the RCWA. Time saving may seem modest and certainly it is if one single point is to be calculated. However, the situation becomes completely different if the code has to be run many times when it is used to build libraries for solving inverse problems.

The other important case of interest is that of crossed gratings where the size of the matrices is squared compared to that of one-dimensional problems. Therefore, within the framework of the RCWA, avoiding solving eigenvalues problems is a significant improvement. Furthermore there is no theoretical difficulty in extending the perturbation method to this latter case. Hence, numerically speaking, we believe that the perturbation method applied to crossed gratings will provide high acceleration ratios.

## Acknowledgments

This work was supported by ANR project LORS.

## References and links

**1. **S. T. Peng, T. Tamir, and H. L. Bertoni, “Theory of periodic dielectric waveguides,” IEEE Trans. Microw. Theory Tech. **MTT-23**, 123–133 (1975). [CrossRef]

**2. **M. G. Moharam and T. K. Gaylord, “Diffraction analysis of dielectric surface-relief gratings,” J. Opt. Soc. Am. A **72**, 1385–1392 (1982). [CrossRef]

**3. **I. C. Botten, M. S. Craig, R. C. McPhedran, J. L. Adams, and J. R. Andrewartha, “The dielectric lamellar diffraction grating,” Opt. Acta (Lond.) **28**, 413–428 (1981).

**4. **I. C. Botten, M. S. Craig, R. C. McPhedran, J. L. Adams, and J. R. Andrewartha, “The finitely conducting lamellar diffraction grating,” Opt. Acta (Lond.) **28**, 1087–1102 (1981).

**5. **K. Knop, “Rigorous diffraction theory for transmission phase gratings with deep rectangular grooves,” J. Opt. Soc. Am. **68**, 1206–1210 (1978). [CrossRef]

**6. **Ph. Lalanne and G. M. Morris, “Highly improved convergence of the coupled-wave method for TM polarization,” J. Opt. Soc. Am. A **13**, 779–784 (1996). [CrossRef]

**7. **G. Granet and B. Guizal, “Efficient implementation of the coupled-wave method for metallic lamellar gratings in TM polarization,” J. Opt. Soc. Am. A **13**, 1019–1023 (1996). [CrossRef]

**8. **L. Li, “Use of Fourier series in the analysis of discontinuous periodic structures,” J. Opt. Soc. Am. A **13**, 1870–1876 (1996). [CrossRef]

**9. **L. Li, “Multilayer modal method for diffraction gratings of arbitrary profile, depth, and permittivity,” J. Opt. Soc. Am. A **10**, 2581–2591 (1993). [CrossRef]

**10. **J. Y. Suratteau, M. Cadilhac, and R. Petit, “Sur la détermination numérique des efficacités de certains réseaux diélectriques profonds,” J. Opt. (Paris) **14**, 273–288 (1983).

**11. **D. M. Pai and K. A. Awada, “Analysis of dielectric gratings of arbitrary profiles and thicknesses,” J. Opt. Soc. Am. A **8**, 755–762 (1991). [CrossRef]

**12. **N. Chateau and J.-P. Hugonin, “Algorithm for the rigorous coupled-wave analysis of grating diffraction,” J. Opt. Soc. Am. A **11**, 1321–1331 (1994). [CrossRef]

**13. **N. P. van der Aa and R. M. M. Mattheij, “Computing shape parameter sensitivity of the field of one-dimensional surface-relief gratings by using an analytical approach based on RCWA,” J. Opt. Soc. Am. A **24**, 2692–2700 (2007). [CrossRef]

**14. **E. Popov, M. Nevière, B. Gralak, and G. Tayeb, “Staircase approximation validity for arbitrary-shaped gratings,” J. Opt. Soc. Am. A **19**, 33–42 (2002). [CrossRef]

**15. **J. Chandezon, D. Maystre, and G. Raoult, “A new theoretical method for diffraction gratings and its numerical application,” J. Opt. (Paris) **11**, 235–241 (1980).

**16. **L. Li, J. Chandezon, G. Granet, and J.-P. Plumey, “Rigorous and efficient grating-analysis method made easy for optical engineers,” Appl. Opt. **38**, 304–313 (1999) [CrossRef]

**17. **T. Vallius, “Comparing the Fourier modal method with the C method: analysis of conducting multilevel gratings in TM polarization,” J. Opt. Soc. Am. A **19**, 1555–1562 (2002). [CrossRef]

**18. **G. Granet, J. Chandezon, J.-P. Plumey, and K. Raniriharinosy, “Reformulation of the coordinate transformation method through the concept of adaptive spatial resolution. Application to trapezoidal gratings,” J. Opt. Soc. Am. A **18**, 2102–2108 (2001). [CrossRef]

**19. **L.D. Landau and E.M. Lifschitz, *Quantum Mechanics: Non-Relativistic Theory*, 3rd ed. Ed. Pergamon, New York, (1977).

**20. **C. Cohen-Tannoudji, B. Diu, and F. Laloë, *Mécanique quantique*, T2, Ed. Hermann, Paris, (1994). [PubMed]