Tryfon Antonakakis

Fadi Baïda

Abderrahmane Belkhir

Kirill Cherednichenko

Shane Cooper

Richard Craster

Guillaume Demesy

John Desanto

Gérard Granet

Editor Evgeny Popov

Boris Gralak

Sébastien Guenneau

Daniel Maystre
email: daniel.maystre@fresnel.fr

André Nicolet

Brian Stout

Gérard Tayeb

Fréderic Zolla

Benjamin Vial

GratinGs: theory and numeric applications

Introduction

Integral methods for scattering problems represent a class of mathematical methods based on integral equations. In this chapter, all the integral equations are deduced from boundary value problems of scattering and are classified as Fredholm integral equations. They can be written in the form: c u(r) v(r) W(r, r )u(r )dr ,

′ ′ ′ = + ∫        (4.1)
in which  may represent for example the surface of a three-dimension diffracting object or the cross-section boundary of a two-dimension (cylindrical) object, v and u are continuous functions in  . The kernel W of the equation is also continuous in  , and c 0 or 1 = according to whether the equation is of the first or second kind. The mathematical problem is to determine u if v and W are known [START_REF] Courant | Methods of Mathematical Physics[END_REF]. This kind of method is widely employed in many domains of physics [START_REF] Courant | Methods of Mathematical Physics[END_REF][START_REF] Morse | Methods of Theoretical Physics[END_REF], where it is usual to extend it to cases in which u, v, and W are only peasewise continuous or can even be singular. A typical example of use in electromagnetism can be found by applying the second theorem of Green for expressing the field in a given volume in terms of its values and of the values of its normal derivative on the surrounding surface. In that case, the kernels of the integral equations include combinations of Green's functions and of their normal derivatives on the boundaries of the scattering objects.

The theory of integral equations can be described in a rigorous, elegant and concise way using distribution theory [START_REF] Schwartz | Méthodes Mathématiques pour les Sciences Physiques[END_REF][4][START_REF] Schwartz | Mathematics for the physical sciences[END_REF] that extends derivatives and other differential operators to discontinuous functions (a famous example is the so-called Dirac function that, for the mathematician, should not be called function). The interested reader can find a detailed presentation of rigorous use of the distribution theory in the electromagnetic theory of gratings in [START_REF] Petit | Electromagnetic Theory of Gratings[END_REF][START_REF] Maystre | Thèse d'Etat[END_REF].

The first application of the integral method in grating theory was proposed almost simultaneously for the case of perfectly conducting gratings by Petit and Cadilhac [START_REF] Petit | [END_REF], Wirgin [9], and Uretski [10]. The first numerical implementation was reported by Petit [11,12] for TE polarization (called also P, or // E , or s polarization). The first extension of this integral equation to the other polarization (TM or S or p) led to a non-integrable kernel. This problem was solved by Pavageau, who proposed for both polarizations new integral equations having continuous and bounded kernels [13].

Soon after, Wirgin [14], Neureuther and Zaki [15], and Van den Berg [16] gave formulations of the integral method applied to gratings made of metals with finite conductivity or dielectrics. The approach was based on the resolution of two coupled integral equations containing two unknown functions. Problems of limited memory storage and timecomputation on computers in the late '60s restricted the numerical implementation of this theory to dielectric gratings only, for which very rare numerical results were published. This restriction was not considered as dramatic by grating specialists at that time. Indeed, it was generally considered that in the visible and infrared regions in which the reflectivity of usual metals exceeds 90%, the model of perfectly conducting grating is accurate. However, development of space optics and astronomy at that time required a precise treatement of the problem of metallic gratings used in the ultraviolet, a domain where the metal reflectivity starts to drop down as approaching electron plasma frequencies. This need required a new approach proposed by Maystre in 1972 [17] by using a single integral equation for a single unknown function. Solving the difficulties in the summation of the series in the kernels and in their integration, the integral method in this formulation was the first one to result in a computer code that was able to correctly model diffraction gratings behaviour over the entire spectrum for almost all commercial gratings profiles [18]. One of the most important conclusions for the practical applications and grating manufacturers was the definite demonstration of the inadequacy of the model of perfectly conducting grating in the nearinfrared, visible and ultraviolet regions [19].

However, the method was unable to treat some kinds of gratings, for example gratings having large periods and steep facets (echelle gratings, for example) or gratings with profile that cannot be represented by Fourier series (rod gratings, cavity gratings, etc.). A further development of this approach was proposed by Maystre [20]. It was numerically implemented in the early '90s in the code 'Grating 2000' by the author, which is used in many academic and industrial centers in the world. Other development was required for other exotic cases that started to find applications and thus required theoretical support for modeling. This development covered conical mountings and specially gratings with dielectric multilayer coatings [21,22], buried gratings and bimetallic gratings [23][24][25][26].

It must be emphasized that in this chapter, the authors use, without complete demonstrations, some analytic properties of gratings demonstrated in chapter 2. Thus, it is recommended to read this chapter before the present one.

First, we will deal with the most frequent problem: the bare metallic or dielectric grating. Then, extensions will be given to other kinds of gratings like perfectly conducting gratings, dielectric coated gratings or gratings in conical diffraction.

4.2.

The integral method applied to a bare, metallic or dielectric grating.

The physical model

The grating surface S of period d separates a region V + with real relative electric permittivity and magnetic permeability + ε and + μ respectively and a region V -with real or complex relative electric permittivity and magnetic permeability - ε and - μ (figure 4.1). The indices + n and n -of these media are given by + + + n ε μ = and n ε μ ---= . We consider the classical diffraction case with incident wavevector i k  lying in the xz plane, i.e. the plane perpendicular to the grooves. The incidence angle i θ is measured in the counterclockwise sense from the z axis and 2 λ= k π denotes the wavelength of light in vacuum. The ordinate of the top of the profile is denoted by 0 z and unit normal S N  is oriented towards V + . We denote by s the curvilinear abscissa on S , with origin being located at the origin of the Cartesian coordinates, and d s denoting the curvilinear abscissa of the point of S of abscissa x d = . In this chapter, we use the complex notation with a time-dependence in ( ) exp iωt -. Let F denote the y-component of the electromagnetic field. In TE polarization, it stays for the ycomponent E y of the electric field, and in TM case for the y-component H y of the magnetic field.

The boundary value problem.

It is shown in chapter 2 that in that case of classical diffraction, the total field

T+ + T T F in V , F = F in V - -     
is invariant along the y axis and that it is pseudo-periodic in x:

( ) ( ) ( ) T T 0 F x+d,z F x,z exp iα d , = (4.2)
with: 0 2 k n sin , k .

+ π α = θ = λ i (4.3)
Moreover, the scattered field is defined by:

T i T F F F in V , F F F in V , + + + - - -  = -  =  =   (4.4)
with the incident field i F in V + being given by: i 0 i x ikn cos z i F (x, z) e . The interest of the notion of scattered field is that it satisfies a radiation condition (also called Sommerfeld condition, or outgoing wave condition, see chapter 2) for z → ±∞ . The radiation condition states that the scattered field at infinity must remain bounded and must propagate upward in + V and dawnward in V -. The scattered field also satisfies Helmholtz equations:

( )

2 2 2 F k n F 0 in V . ± ± ± ± ∇ + = (4.6)
The invariance along the y axis allows one to reduce the scattering problem to a twodimension problem while the pseudo-periodicity restricts the study of the field to a single period of the grating. Consequently, it can be considered that V ± are no more volumes but surfaces extending on a single period of the grating.

In order to periodize the scattered field, we introduce the function

+ + U in V , U= U in V , - -      : ( ) ( ) 0 i x U x
, z e F x, z .

-α = (4.7)

We denote by ( ) The Helmholtz equations and the radiation condition are not sufficient to define the boundary value problem satisfied by the scattered field F . A third kind of condition must be added: the boundary conditions of the electromagnetic field components across S. The tangential components of the electric and magnetic fields are continuous across this interface, as far as the permittivities and permeabilities of the two media take finite values. For both polarizations, this property yields: [ ] i i 0 F x(s), z(s) exp( i x(s)) exp ikn cos z(s) .

ψ
+   ψ = -α = - θ   i (4.9)
Using Maxwell equations, the continuity of the tangential component of the magnetic field (for TE polarization) and that of the electric field (for TM polarization) leads to the following relation: i q (s) (s) q (s),

+ + --   φ + φ = φ   (4.10) with ( ) ( ) ( ) i i 0 S F exp( i x)
N dz s dx s ikn sin cos exp ikn cos z(s) , ds ds

+ + ∂ φ = -α = ∂   =- θ + θ - θ     i i i (4.11)
and 1 , for TE polarization, q 1 , for TM polarization.

± ± ±   µ  =    ε  (4.12)

Integral equation

The theoretical basis of the integral method lies on a general property of the electromagnetic field: the field inside a given surface of the xz plane can be expressed from the values of the field and of its normal derivatives on the curve surrounding the surface, according to the second Green's theorem. The value of U ± at a point of V ± of coordinates x and z be deduced from its values on S using equation (4.139) of appendix 4.A:

d s s' 0 U (x, z) (x, z,s ') (s ') (x, z,s ') (s ') ds ', ± ± ± ± ± =   = ± φ + ψ   ∫ G N (4.13) with ( ) ( ) { } m m m 1 1 (x, z,s ') exp imK x x ' s ' i z z ' s ' , 2id ∞ ± ± ± =-∞ = - + γ -     γ ∑ G (4.14) ( ) ( ) ( ) { } m m m m 1 dx(s ') dz '(s ') (x, z,s ') sgn z z ' s ' 2d ds ' ds ' exp imK x x ' s ' i z z ' s ' , ∞ ± ± =-∞ ±   α   = - - ×       γ     × - + γ -     ∑ N (4.15)
where:

m 0 mK, α = α + (4.16) 2 K , d π = (4.17) ( ) 2 2 m m kn , ± ± γ = -α (4.18)
with s' being the curvilinear abscissa on a point of S with coordinates x '(s ') and z'(s') .

According to section 4.A.

+ + + + + = ψ   φ + ψ - =   ∫    G N , (4.19) d s s' 0 (s ') (s,s ') (s ') (s,s ') (s ') ds ' 0 2 - - - - - = ψ   φ + ψ + =   ∫    G N , (4.20)
with:

( ) ( ) ( ) ( ) ( ) m m m 1 1 (s,s ') exp imK x s x ' s ' i z s z ' s ' 2id ∞ ± ± ± =-∞   = - + γ -   γ ∑ G , (4.21) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) m m m m 1 dx ' dz ' (s,s ') sgn z s z ' s ' 2d ds ' ds ' exp imK x s x ' s ' i z s z ' s ' . ∞ ± ± =-∞ ±   α = - - ×   γ       × - + γ -   ∑ N (4.22)
Introducing in eq. (4.20) the values of i (s) (s) (s)

- + ψ = ψ + ψ and i q (s) (s) (s) q + - + -   φ = φ + φ  
given by the continuity conditions on the grating profile (eqs (4.8), (4.9), (4.10), (4.11) and (4.12)) yields a second integral equations with two unknown functions + ψ and + φ :

d s s' 0 i i i (s,s ') (s,s ') ds ' 0. 2 q (s ') (s ') (s ') (s ') q (s ') (s ') - - = + + + - +     φ + φ ψ + ψ       +      ψ + ψ    + = ∫ G N (4.
+ + - -γ + + + =     -φ α = + - ψ       γ γ       ∫ (4.24) ( ) ( ) d m s imKx s i z s m m m m s 0 i (s) 1 dx(s) dz(s) t e (s) ds, 2d ds ds - - - + γ - - + =     φ α = + + ψ       γ γ       ∫ (4.25)
with 0 z being the ordinate of the top of the grating profile. For non-evenescent reflected orders, diffraction efficiencies m ρ can be obtained using eq. (4.187):

2 m m m 0 r . + + γ ρ = γ (4.26)
For gratings made of a lossless dielectric material in V -, transmitted efficiencies can be defined as well:

2 m m m 0 q t . q + - + + γ τ = γ (4.27)
In that case the energy balance (see chapter 2) can be expressed by: m P m P m m 1,

+ - ∈ ∈ ρ + τ = ∑ ∑ (4.28)
with P + and P -denoting respectively the set of non-evanescent reflected and transmitted orders. The numerical implementation of the integral equations will be described in section 6.

In contrast with the two coupled equations obtained in this section, the integral equation obtained by Maystre is unique. This feature requires the definition of a single and well adapted unknown function. The mathematical definition of this function needs the use of tools of applied mathematics described in appendix A. Appendix B contains a mathematical description of this mathematical function and of the integral equation. Here, we give a heuristic description of this function for TE polarization. First, we replace the material in V - by the same material as in V + , the entire space being thus homogeneous. It can be shown that it exists one (and only one) distribution of surface current density Φ parallel to the y axis, placed on S (this surface separates now two identical media), which generates in V + a field equal to the actual diffracted field in the physical problem. Intuitively, it is easy, from this surface current density, to express in an integral form the actual scattered field in V + , this current distribution being nothing else than a set of current lines placed in a homogeneous medium. From the expression of the scattered field in V + , simple mathematical calculations allow one to deduce the scattered field and its normal derivative above S, thus

φ thus, using the second Green theorem, of the actual physical field below S. At that point it has been shown that the four unknown functions contained in the classical theory previously described in this section can be derived from a single one and that, in some way, there is a redundancy in the use of multiple unknown functions. It is easy to understand that this single unknown function can be calculated from a single integral equation. This equation can be obtained for example by writing the continuity on S of the integral expressions of the field in V + and V -.

This method was the first one to show that, in contrast with the second Green theorem, it exists a formula that allows one to express the field inside a given domain from a single function defined on its boundary. This function is neither the field nor its normal derivative, but both can be deduced from it through simple integrals. These integrals automatically satisfy the compatibility condition.

The bare, perfectly conducting grating

Perfectly conducting gratings were historically the first gratings to be modeled using rigorous electromagnetic theories. They represent accurate models for metallic gratings working in far infrared and microwaves regions. The pioneering works appear in the '60s [START_REF] Petit | [END_REF] and were followed by many papers. Various formulations of the integral method have been published. They differ either in the form of the integral equation ot in the numerical implementation. A review of this matter may be found in [27,[START_REF] Petit | Electromagnetic Theory of Gratings[END_REF][START_REF] Maystre | Rigorous vector theories of diffraction gratings[END_REF].

Perfectly conducting gratings in TE polarization

Two different approaches will be presented in this section. The first one, published by R. Petit and M. Cadilhac in [START_REF] Petit | [END_REF], leads to a Fredholm integral equation of the fist kind with a singular kernel. A version leading to a Fredholm integral equation of the second kind with a nonsingular kernel was proposed by Pavageau et al. [13] using the ideas of Maue [29].

The total field in V + is pseudo-periodic, it satisfies the Helmholtz equation:

( ) 2 2 T 2 T F k n F 0 in V , + + + + ∇ + = (4.29)
and a radiation condition at infinity. Thus we can apply the generalized compatibility condition of section 4.A.5:

( ) d s i s' 0 (s) (s,s ') (s ') (s,s ') (s ') ds ' s , 2 + + + + + = Ψ   Φ + Ψ + ψ =   ∫ G N (4.
m m m 1 1 (s,s ') exp imK x s x ' s ' i z s z ' s ' . 2id ∞ + + + =-∞   = - + γ -   γ ∑ G (4.32) ( ) ( ) ( ) ( ) ( )
This is a Fredholm integral equation of the first kind, with a singular kernel. The amplitudes of the reflected waves are deduced from eq. (4.186):

( ) ( ) ( ) d s m m m s 0 1 r
exp imKx s i z s s ds 2id

+ + + =   = - -γ Φ   γ ∫ . ( 4
+ -+ γ - + =-∞ +∞   α = - -   γ     ∑ K . ( 4
+ + + + + = Ψ   Φ + Ψ + ψ =   ∫ G N (4.37)
In that case too, the tangential component of the electric field vanishes on the profile. It is shown in chapter 2 from Maxwell equations that this condition entails that the normal derivative of the total field vanishes on S, thus :

0 + Φ =, (4.38) so that ( ) d s i s' 0 (s)
∞ + ± =-∞ +   α = - - ×   γ       × - + γ -   ∑ N (4.40)
The Fredholm integral equation of the second kind with a regular continuous kernel is very close to that obtained for TE polarization (eq. (4.35)). The amplitudes of the reflected waves are deduced from eq. (4.186):

( ) ( ) ( ) d s m m m m s 0 1 dx(s) dz(s) r
exp imKx s i z s s ds 2d ds ds

+ + =   α   = - -γ - Ψ       γ   ∫ , (4.41)
As for TE polarization, the efficiencies

Multiprofile gratings

The use of dielectric coatings has many applications even for diffraction gratings use. For example, metallic gratings are covered by a thin layer of dielectric material in order to avoid oxidation of the metal. Dielectric gratings can require an antireflection coating consisting of a thin layer or a stack of layers. Conversely, a stack of layers is used to increase the metal reflectivity or even to replace it, in order to reduce the absorption of light beams for high power laser applications. Up to our knowledge, the first numerical results in the study of coated gratings was made by Van de Berg [16] for a single-layer perfectly conducting grating with the layer filling up the space between the grating surface and a plane surface. Although at that time the interest in such geometry remained mostly academic, further development of technology made it possible to fabricate such layer by dielectric coating and polishing. Another important application comes from the process of replication of dielectric gratings using an epoxy layer to transfer the replica to a plane surface of the substrate, or to have epoxy as grating layer itself.

After this initial work, two integral methods were proposed. The first one [20] is theoretically able to deal with an arbitrary multilayer grating without limitations concerning the shape of the profile or the conductivity of the layers. The second method [21] can deal with a multilayer grating without interpenetration of the profiles.Botten has solved the problem with a single-profile grating that has a stack of plane layers below and, eventually, above it [22,23], by introducing a new form of Green's function, adapted to the multilayer system, which leeds to a single integral equation.

In what follows, we will at first describe the method for a single interface inside a stack, when the layer is relatively thin so that the upper and the lower interface interpenetrate. It is important to distinguish the two cases, with and without interpenetration, because in the latter case, it is possible to define a plane layer in between that does not cross the upper or the lower interface. This possibility enables one to use the plane-wave Rayleigh expansion of the electromagnetic field between the interfaces, whereas in the former case it is necessary to write and to solve a system of coupled integral equations that link the field components on the top and bottom of each layer. Let us consider the case of a multilayer grating having profiles that interpenetrate. In other words, it is impossible to introduce inside the layer a plane surface that does not cross one of the profiles. Then it is impssible to use the plane-wave expansion between the interfaces and we are led to solve integral equations that are coupled on the two interfaces of each layer.

Thin-layer gratings

We introduce in figure 4.2 a grating made of M materials (numbered from 0 to M, separated by M-1 profiles (numbered from 1 to M). We introduce the following functions:

V j-1 V j V j+1 S j S j+1
D. Maystre and E. Popov: Integral Method for Gratings 4.11

0 0 i x i x i j j,M j j F e F e in V U . 0 elsewhere. -α -α  - δ  =    (4.43)
The function U j inside each layer can be expressed from the fields and normal derivatives on the lower and upper interface, for j=1, M-1: j j j 1 j 1 j j j j j j j j j j j S S j 1 j 1 j 1 j 1 j 1 j 1 j 1 j 1 j 1 j 1 S S U (x, y) (x, y,s ) (s )ds ' (x, y,s ) (s )ds (x, y,s ) (s )ds (x, y,s ) (s )ds .

+ + + + + + - - - - + + + + + + + + + + ′ ′ ′ ′ ′ = φ + ψ ′ ′ ′ ′ ′ ′ - φ - ψ ∫ ∫ ∫ ∫ G N G N (4.44)
The expression being limited to the second one if j = 0 and to the first one for j = M. The functions derived from the Green functions in the various materials depend on the interface number:

( ) ( ) { } j j j,m j m j,m 1 1 (x, y,s ') exp imK x x ' s i z z ' s 2id ∞ ± ± ± =-∞   ′ ′ = - + γ -   γ ∑ G , (4.45) ( ) ( ) ( ) {
} j m j j m j,m j j,m j dz '(s ) 1 dx '(s ) (x, y,s ) sgn z z ' s 2d ds ' ds ' exp imK x x ' s i z z' s , ∞ ± ± =-∞ ±   ′ ′ α     ′ ′ = - - ×     γ       ′ ′ × - + γ -   ∑ N (4.46) with:
( )

2 2 j,m j 1,m j m kn + - + γ = γ = -α . (4.47)
It can be shown it the same manner as in eq. (4.146) that a compatibility condition on the j th interfacewritten in a matrix form is given by: j j j j j j,j 1 j 1 j,j 1 j 1 2

+ + + + + - - - - + + + + ψ = φ + ψ - φ - ψ G N G N . (4.48)
Another compatibility equation is obtained on the (j+1) th interface: j 1 j 1,j j j 1,j j j 1 j 1 j 1 j 1 2

- + + + - + -- - - + + + + + + ψ = φ ψ - φ - ψ + G N G N . ( 4.49)
In eqs. (4.48) and (4.49), we use the double-index Greens functions that are derived on two consecutive profiles: j,j 1 j j j 1 (s ,s )

- - + + ′ ≡ G G , j 1,j j j 1 j
(s ,s )

+ + + + ′ ≡ G G
, and similarly for N .

The computability equation becomes, for the upper and lower media:

M M M M M 2 + + + + + ψ = φ + ψ G N , (4.50) 1 1 1 1 1 2 - -- -- ψ = -φ - ψ G N . ( 4.51)
By combining eqs (4.48) with eq. (4.49), we obtain the link between the unknows on the upper and on the lower interface of the j th layer, for j=1, M-1: j, j 1 j, j 1 j j j 1 j j 1 j 1 j 1 j j 1, j j 1, j 2 2

- - + + - + + + + - - - + + + + + + + +         - ψ ψ                 + φ φ               = N G N G N G N G . (4.52)
This equation gives the transmission operator of the unknown amplitudes across the j th layer, i.e. from the j th to the (j+1) st interface, for j=1,M-1:

1 j, j 1 j, j 1 j j j 1, j j 1 j 1 j 1, j j 1, j 2 T 2 - - - + + + + + - - + + + + + +     -         +           = N G N G N G N G . (4.53)
The transmission operator includes an inverse operator. Numerically, this inversion leads to the inversion of a matrix, as we will see in section 4.6.

The transmission matrix across the j th interface for j = 1, …, M-1 is obtained using the continuity of the tangential and normal field components, as given by eqs.(4.8) and (4.10): j j j j j j j j 0 q T , T 0 q

+ - - +- +- + - +       ψ ψ       = =         φ φ           , j M ≠ . (4.54)
The advantage of this presentation is that there are no exponentially growing terms in the transmission matrices, since all the components of the two variable functions contain only scattered propagating or decreasing evanescent waves, in contrast with the other methods (differential, Fourier modal, Rayleigh, etc.). However, this formulation requires calculating the cross-layer functions between the interfaces, which leads to computation times equal to those of single-interface functions. As a consequence, the total computation time is almost multiplied by a factor 2 with respect to the case where cross-layer kernels can be avoided, as discussed in the next section.

Finally, it can be deduced from eqs. (4.53) and (4.54

): i M 1 M M,M 1 2 2,1 1 i M 1
T T ... T T T .

+ - +- +- +- - + -       ψ ψ ψ       + =       φ φ φ       (4.55)
Finally, eqs. (4.55), (4.50) and (4.51) form an operator system of 4 equations with 4 unknown functions, which can be solved on a computer after representing each operator by a matrix, as described in section 6.

Profiles without interpenetration

This case is simpler than the situation in sec.4.1, because it is possible to use the plane-wave expansion between the grating profiles and thus to decouple the integral presentation used in eq.(4.44). The idea is illustrated in figure 4.3.

If it is possible to introduce a plane layer between the profiles, the plane wave expansion is valid inside this layer. The advantage is that the plane waves (propagating and evanescent) that represent each diffraction order m can be easily separated into to sets: (i) upgoing waves having amplitudes of the y-component of the field equal to j,m r that are generated by the lower surface S j , (ii) downgoing waves with amplitudes j,m t generated by the upper grating surface S j+1 .

+ + + + = ψ + φ   =   ∫ , (4.56) with ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )
j,m j j j,m j j,m j j m j,m j j j,m j j j j,m

+ + + + + +   = - -γ   γ   α     = - - -γ     γ   (4.57)
and with j s denoting the curvilinear abscissa on the j th profile, j,d s being the curvilinear abscissa of the point of S j of abscissa x = d.

We then can directly use eq.(4.48) to express the field on the interface z = j 1 z(s ) + :

( )

j 1 j 1 j 1 j,m j 1 j,m j 1 j 1 j 1 j 1 m s exp imKx(s ) i z(s ) r 2 - + + + + -- - - + + + + + + ψ   = + γ - φ - ψ   ∑ G N . ( 4

.58)

Let us consider the sum in eq. ( 4.58):

j 1 j 1 j,m j 1 j,m m exp imKx(s ) i z(s ) r + + + +   ζ = + γ   ∑ . (4.59)
Inserting in this equation the value of j,m r given by eq. (4.56) yields:

x z j z + j z - V j-1 V j V j+1 S j S j+1 ( ) ( ) ( ) ( ) j,d j j 1 j 1 s j 1 j,m j 1 j,m j j j j,m j j j j m s 0 (s ) exp imKx(s ) i z(s ) N s s G s s ds , + + + + + + + + + = ζ =     + γ ψ + φ     ∑ ∫ (4.60)
which can be written in operator form: j 1 j,j 1 j j,j 1 j ,

+ + + + + + + ζ = ψ + φ N G (4.61)
where the operators j, j 1 + + N and j, j 1

+ +
G are obtained from a summation in m and an integral in j s . Thus we can write eq. ( 4.58) in the operator form: j 1 j,j 1 j j,j 1 j j 1 j 1 j 1 j 1 . 2

- + + + + + -- - - + + + + + + ψ = ψ + φ - φ - ψ G N N G (4.62)
The second relation comes from eq. (4.49) by using the amplitudes of the diffraction orders coming down from the upper interfaces and valid below z = j 1 z - + . Following the same lines as in the previous paragraph yields:

j j 1,j j 1 j 1,j j 1 j j j j 2 - - - - - + + + + + + + + ψ = ψ + φ φ ψ N G + G + N . ( 4.63)
The transmission matrix between the j th and the j+1 st interface takes a form similar to the case of interpenetrating layers, eq.( 4.53). However, the difference is essential, because each series used in (4.63) is evaluated on a single interface:

1 j 1, j j 1, j j j j 1, j j 1 j 1 j, j 1 j, j 1 2 T 2 - - - + + + + + - - + + + + + +     -         +           N G = N G N G N G . (4.64)
The second difference is that the exponential terms are explicitly given in the j, j 1

+ + N , j 1, j - + N , j, j 1 + + G and j 1, j - + G through the functions j 1 j,m j 1 exp imKx(s ) i z(s ) + + +   + γ   and j j,m j exp imKx(s ) i z(s ) +   -γ  
in such a way that it can be extracted from each of these operators a part containing all the growing and decreasing exponential terms, which allow a much better stability of the numerical implementation through adequate treatments, for example the S-matrix algorithm described in appendix A and B of Chapter 7.

Gratings in conical mounting

When classical gratings with one-dimensional periodicity are used with incidence plane perpendicular to the grooves, the diffraction orders lie in the same plane. Off-plane incidence brings the diffraction orders out of the plane and their directions lie on a cone , which explains the term of conical diffraction. One of the first experimental works can be found in [24,25]. The use of conical mount is typical for concave gratings and in some spectrographs aiming to separate off-plane the incident and the diffracted directions.

The first theoretical studies were made in 1971 using the integral [26] and the differential [30] methods. An interesting theoretical development came in 1972 when Maystre and Petit demonstrated that under special conditions, mechanically ruled perfectly conducting gratings can have very high and constant efficiency over a large spectral domain [31]. They also established the theorem of invariance [32] that gives an expression of the diffraction efficiency of a perfectly conducting gratings in conical mounting, expressed as a linear combination of efficiencies in an in-plane (classical) mount for the two fundamental polarizations. Since the theorem is not valid for finitely conducting metals, later development of the integral method allowed studies of diffraction gratings behaviour in conical mount when working in the UV and visible [33,34]. An interested reader can find the demonstration of the invariance theorem in [START_REF] Maystre | Thèse d'Etat[END_REF]33]. (ii) fictitious equivalent classical case:

The wavelength is increased to / cos , ′ λ =λ ϕ the angle ′ ϕ is now equal to 0 (in-plane incidence), the angle ′ θ = θ i i , then efficiencies in TE and TM polarizations are denoted by ( )

ΤΕ i m ρ , ′ θ λ and ( ) ΤM i m ρ
, ′ θ λ , respectively. It can be noticed that the projection i xz k  of the wavevector of the incident wave on the xz plane in conical mount identifies with the wavevector of the incident waves in the fictitious equivalent case.

The invariance theorem states that:

( ) ( ) ( ) ( ) ( ) 2 2 i ΤΕ i ΤM i m m m , , , cos ρ , sin ρ , ′ ′ ρ θ ϕ δ λ = δ θ λ + δ θ λ . (4.65)
It is to be noticed that, as for the incident wavevector, the projections of the wavevectors of the scattered waves on the xz plane in conical mount are identical to the wavevectors of the scattered waves in the fictitious equivalent case.

Numerical tools for an efficient numerical implementation

Integration schemes for the integral equation

All the integral equations in this chapter link the value of an unknown function u(s) at a given point of S to its value on its entire domain of definition:

d s 0 c u(s) v(s) W(s,s ')u(s ')ds ' = + ∫ , ( 4.66)
where all functions are periodical, with v and W being known functions, W being possibly singular but integrable. The constant c takes values 0 or 1 according to whether the integral equation is of the first or second kind.

There are many different ways to solve such an equation for the grating problem. Pavageau et al. proposed an iterative method [13] that does not require any matrix inversion, like the well known Born method for scattering problems. Unfortunately, it may diverge [35,36].

A well known general methods is based on the periodicity of all functions of the equation, which allows a projection of these functions and of the equation on the Fourier space:

m s s d m 2 u(s) u exp(imK s), K , s π = = ∑ (4.67)
and similar expressions for v and W. The integral equation is transformed into a linear system of algebraic equations:

( ) nm nm m n m W c u v , n -δ = ∀ ∑ , ( 4.68)
which can be solved numerically after truncation. However, this approach requires computing a double Fourier decomposition:

d d s s nm s s 0 0 W W(s,s ')u(s ') exp( inK s imK s ')ds ds ' = - - ∫ ∫ . (4.69)
The method has been applied to gratings with profiles consisting of few straight segments, because in that case the double Fourier integral can be calculated in closed form. It is so for triangular profiles [11] or trapezoidal gratings [37].

The most widespread method is the so-called point-matching (or discretization) method. Instead ot using discrete Fourier components, the unknown function is discretized on the grating profile and represented by its vales j j u u(s ) = inside the interval of integration. Ther integration process leads toan equation quite similar to eq.(4.68) ( ) (multiplied by the weight of integration) and can require much more complicated treatments, specially if W is singular. The sophistication of this treatment is one of the decisive keys for the precision of the solution of the integral solution. The second key is the analytical study of the kernel, which is described in the next two sections. When W is regular, continuous, with a continuous first derivative, the rectangular rule of integration is quite precise since the function to be integrated is periodic and it can be noticed that the rectangular and trapezoidal rules are completely equivalent in that case. However, the derivative of W(s,s') is generally discontinuous when s = s' and a trapezoidal rule or higher order treatment provide a better precision [38,39]. In what follows we assume that u(s) is a continuous function. This is obviously the case when the grating profile has no edges. Several more detailed arguments in favour of the rectangular rule can be found in [START_REF] Petit | Electromagnetic Theory of Gratings[END_REF][START_REF] Maystre | Thèse d'Etat[END_REF]. Let us shortly repeat one of them. We consider an integral of a periodical continuous function a(s): Simple calculations show that when using the rectangular rule with P discretization points, the integration error is proportional to the P th Fourier coefficient of a(x). Since it is continuous, the Fourier series converges like 1/P 2 at least.

jp jp p j p W c u v , j [1, P] -δ = ∀ ∈ ∑ , ( 4
The implementation of the rectangular rule is very simple. We define the values of jp W in the following manner:

d j,p d d s j p W W s , s , j, p 1,..., P 1 P P P   = = -     . (4.73)
Using this result, the product of the unknown functions u(s) with the non-singular parts of the kernels can be integrated in the form of a simple matrix product: In the case of multilayer gratings with profiles interpenetration, there are two main situations that can complicate the numerical evaluation of the cross-layer functions which link the fields and normal derivatives on both sides of a layer: (i)

The wavelength λ is much larger than the layer thickness t. In that case the profiles are located too close to each other, compared to λ, so that the functions G and N become large in modulus for j 1 j s s + → . This behaviour has the same origin as the singularities of the kernels G and N for a bare grating , which will be discussed in section 6.3 and can be eliminated in a similar manner. Numerical results show no problems as long as layer thicknesses exceed λ/20.

(ii)

The period d is much greater than the layer thickness t. If the wavelength is not too much larger than t, then the kernels have no singularities, but their moduli present peaks when the distance j 1,p j,q P P + between two points located on the two different profiles is small with respect to the discretization distance between the points located on the same profile. The width of these maxima is of the order of magnitude of the layer thickness, thus the correct implementation the trapezoidal integration rule requires the distance between two consecutive points of the profile discretization j,p j,p 1 s s - ∆ = to be less than the width of the maxima.

As a rule of thumb, if

p d N ∆ ≈
, where N p is the number of integration points, then its lower limit is determined by the relation:

p,min d N t ∝ .
(4.75) Thus, for echelles, a values of d of about 10 µm and t of 20 nm requires the number of integration points to exceed 500. Note that N p directly determines the number of unknown values of φ j and ψ j and thus the size of the matrices to be used. Practically, it is not worth nowadays for N p to exceed 10 000, because of the computation time, memory requirements, round-off errors and limited digits. As a consequence, it is necessary to find another way of integration instead of the trapezoidal rule.

There is another problem that can come from the matrix inversion in the construction of the transmission matrix between the layers, eq.(4.53). Contrary to the transmission matrix that contains growing and decreasing exponential terms in the plane wave expansion, (thus requires some type of recursive algorithm to contain the contamination of the growing exponentials, S-matrix algorithm, for example, see appendix 4.C), the distance between the profiles in the z-direction that appears in the kernels in eqs. (4.45) and (4.46) restricts the terms to only propagation or evanescently decreasing ones. However, the matrix inversion of these terms that is required in eq.(4.53) can create exponentially growing terms. Two situations can appear: 1. The matrix inversion in eq.(4.53) can be done without numerical problems. This happens when the layer thickness is not quite large. In that case it is possible to progress upwards in the stack of layers by following the S-matrix algorithm. 2. The matrix inversion does not work. This could happen if the distance between two consecutive interface is large. Two different geometries can be concerned:

(i)
there is no interpenetration of these two profiles. In that case one can easily apply the technique described in the next section. (ii) there is interpenetration of two very deep interfaces. It is possible to use directly eq.(4.52) in the S-matrix algorithm without inverting the matrix to calculate the entire Tmatrix in eq. (4.53). The formulation of the S-matrix algorithm to an equation having the form given in (4.52) is quite similar to the classical aplication, but it requires one additional iteration step. The advantage is that it avoids the inversion of small terms that can lead to singular matrices. This special technique is given in Appendix 4.C and is not quite popular, but can be used in other methods that apply for multilayer stack, for example, in the coordinate transformation method.

Summation of the kernels

There are several problems in the calculation of the functions included in eqs (4. 19) and (4.23):

( ) ( ) ( ) ( ) m imK x s x ' s' i z s z' s' m m 1 1 (s,s ') e 2id ± ∞ - + γ -   ±   ± =-∞ = γ ∑ G , (4.76) ( ) ( ) ( ) ( ) ( ) ( ) ( ) m imK x s x ' s' i z s z' s' m m m 1 dx ' dz ' (s,s ') sgn z s z ' s ' e , 2d ds ' ds ' ± ∞ - + γ -   ±   ± =-∞   α = - -   γ     ∑ N (4.77) ( ) ( ) ( ) ( ) ( ) m imK x s x ' s' i z s z' s' m m , m 1 dx dz (s,s ') sgn z z ' e , 2d
ds ds

± - + γ -     + =-∞ +∞   α = - -   γ     ∑ K (4.78) with: m 0 mK α = α + , (4.79) ( ) 2 2 m m kn ± ± γ = -α . (4.80)
For the sake of simplicity, we assume here that n 1 ± = and we cancel the superscript ± in

ik i i mK , 2 1 1 k 1 . i i m K 2i →∞ →∞ →∞ α → γ ≈ α - → α + α ≈ + → γ α α (4.81) We consider the function [ ] { } m m m 1 1 (s,s ) exp imK x(s) x(s ) i z(s) z(s ) 2id ∞ =-∞ ′ ′ ′ = - + γ - γ ∑ G . (4.82)
At point s s′ = , it is obvious that the sum does not converge since the terms decrease in 1 m .

Of course, a very slow convergence can be expected when the two points are close to each other. Different techniques have been proposed to accelerate the convergence. Neureuther and Zaki [15] have employed a transformation technique based on the use of Mellin transforms.

Other authors have proposed accelerating processes [40][41][42][43]. Here we describe a direct approach [START_REF] Maystre | Thèse d'Etat[END_REF]. Let us determine at first the asymptotic expression of the kernel. If we replace (4.81) into eq.( 4.82), we obtain the asymptotic term in the sum:

[ ] [ ] 0 0 mK x(s) x(s ) i z(s) z(s ) x(s) x(s ) m 1 mK x(s) x(s ) i z(s) z(s ) x(s) x(s ) m 1 1 (s,s ) e e 4 m 1 e e , 4 m -∞ ′ ′ - + -   ′ α -   ∞ =- ∞ ′ ′ - - - -   ′ -α -   = ′ = π - π ∑ ∑ G (4.
0 0 K x(s) x(s ) i z(s) z(s ) x(s) x(s ) K x(s) x(s ) i z(s) z(s ) x(s) x(s ) 1 (s,s ) e log 1 e 4 1 e log 1 e . 4 ′ ′ - - + -   ′ α -   ∞ ′ ′ - - - -   ′ -α -     ′ = -   π     + -   π   G (4.84)
The calculation of the kernel in (4.82) is achieved by subtracting the asymptotic value:

( ) ∞ ∞ = + - G G G G . (4.85)
The first term in the right-and side is explicitely given in (4.84). It is singular for s = s'e. This singularity is integrable and is be treated in the next section.

The term between parenthesis in eq.( 4.85) is obtained by combining the terms in the sums in eqs.(4.82) and (4.83). As far as the second one is the asymptotic value of the former, the series converges, whatever the values of s and s'. Furthermore, it is possible to show that by combining the terms m and -m in the sum, we finally obtain a rapidly converging series, whose terms decrease as m -3 when s s′ = , as shown later in eq.( 4.87) Moreover, in this case this series is continuous and its value is simply given by: ( )

0 m s s m 0 1 1 1 2id 2id 4 m ∞ ′ = ≠   - = + +     γ γ π   ∑ G G . ( 4

.86)

Using the third identity of eq.(4.81), for large values of m, the term in the sum is equal to:

2 2 2 3 3 m m k 4d 1 1 d 1 2id 4 m 8 m + → - → - γ π πλ α . (4.87)
Obviously, the singularity and the slow convergence of G have been carried out by introducing the series ∞ G . Fortunaley, this singularity is logarithmic and thus integrable, as shown in the next section.

Let us now deal with the second function defined in eq. (4.77):

[ ] [ ] m imK x(s) x(s ) i z(s) z(s ) m m m 1 dx(s ) dz(s ) (s,s ) sgn z(s) z(s ) e . 2d
ds ds

∞ ′ ′ - + γ - =-∞   ′ ′ α   ′ ′ = - -   ′ ′ γ     ∑ N ( 4
∞ ∞ ′ ′ ′ -α - - - + - = ′ ′ ′ α - - + - ′ ′ = - ′ ′   ′ - -   ′ ′   ′ ′   ′ + - +   ′ ′   ∑ N 1 -∞ =- ∑ (4.90)
The sums can be evaluated in a closed form:

[ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] 0 0 x(s) x(s ) K x(s) x(s ) imK z(s) z(s ) x(s) x(s ) K x(s) x(s ) imK z(s) z(s )
∞ ′ -α - ′ ′ - - - ′ α - ′ ′ - + - ′ ′ = - ′ ′   ′ + - -   ′ ′   - ′ ′   ′ + - +   ′ ′   - N (4.91)
As noticed for G , the term ∞ -N N must be considered as a series each term of which is obtained by making the difference of the corresponding terms in the sums in eqs.(4.88) and (4.90). This series converges much more rapidly than the series in eq. ( 4.88) and it is continuous at s = s'.

The limit of N when s ' s → can be determined calculating the limits of the two terms ∞ N and ( ) ∞ -N N . After tedious calculations, we can deduce that this limit exists and is given by:

2 2 m 2 2 s' s 0 m m d z dz i 1 1 ds (s,s) lim (s,s ) ds 2 2d 4 dx dz ds ds ∞ → =-∞   α ′ = = - + +     πα γ π       +         ∑ N N . ( 4

.92)

This interesting results established by Pavageau and Bousquet [44] is very important for the numerical applications, because it shows that (s,s) N contains a series that canverges like 1/m 3 (after adding the terms with negative and positive values of m).

Let us notice that the second derivative of the profile function appears in eq.( 4.92) and it clearly requires the continuity of the first derivative, i.e., the absence of edges.

The third kernel (s,s ') K , given by eq. (4.36):

[ ] [ ] m imK x(s) x(s') i z(s) z(s') m m , m 1 dx dz (s,s ') sgn z(s) z(s ') e 2d
ds ds

+ - + γ - + =-∞ +∞   α = - -   γ     ∑ K (4.
∞ → =-∞   α ′ = = - + -     πα γ π       +         ∑ K N . ( 4
.94)

Integration of kernel singularities

Clearly, the asymptotic part of (s,s )

∞ ′ G
in eq.( 4.84) is singular when s s′ → . After some calculations, it can be found that

[ ] [ ] { } ( ) s s 2 2 2 1 (s,s ) ln K x(s) x li (s ) z(s) z(s ) 4 m ′ → ∞ ′ ′ ′ = - + - π G . ( 4
∞ ′ →     ′ -       ′ = π + + +         π             G . (4.96)
The first two terms represent regular parts that can be integrated by using the rectangular rule, as shown later. Unfortunately, s s ln d

′ -

is not periodic and the rectangular rule is very poor when applied to nonperiodic functions. It is possible to overcome this difficulty [START_REF] Maystre | Thèse d'Etat[END_REF] by considering another function defined on (0, d): . The advantage of his new function is that it is continuous except on the singularity, and all its derivatives with respect to s are the same on s′ = 0 and s′ = d and thus we can use the rectangular rule.

s s s s 1 (s,s ) ln ln 1 , s,s (0, d) 2 d d ∞   ′ ′ -  - ′ ′ = + - ∈     π        G , ( 4
We perform the integration of (s,s )

∞ ′ G
by setting:

d d d s s s 0 0 0 (s,s ) (s )ds (s,s ) (s ) (s,s ) (s) ds (s,s ) (s)ds ∞ ∞ ∞ ∞   ′ ′ ′ ′ ′ ′ ′ ′ ′ φ = φ - φ + φ   ∫ ∫ ∫   G G G G . ( 4

.98)

The term in the square bracket is non-singular and can be integrated using the rectangular rule, if we take into account that:

2 2 s s 1 1 dx dz lim (s,s ) (s ) (s,s ) (s) ln(2 ) ln (s) 2 2 ds ds ∞ ∞ ′→             ′ ′ ′ φ - φ = π + + φ           π              G G . ( 4

.99)

The second term in eq.( 4.98) contains the singular part, but it can be integrated analytically:

d d s s 0 0 d (s,s ) (s)ds (s) (s,s )ds (s) ∞ ∞ ′ ′ ′ ′ φ = φ = -φ π ∫ ∫   G G . (4.100)
In conclusion, the integration of the singular kernel is made by introduction at first ,

∞

G

which permits to define a series ∞ -G G that is continuous and rapidly converging hence easily integrable by the use of the rectangular rule. Second, the integration of the term containing ∞ G is performed by defining a new function ∞  G , which has the same singularity as ∞ G , has the property of a periodic function, and can be analytically integrated.

Kernel singularities for highly conducting metals

When the conductivity of a metallic grating tends to infinity, the Green function in the metal tends to a delta function. This property is rather obvious: for very large conductivities, the field generated by a line current (delta function) placed in the metal decreases very rapidly since it is absorbed on very short distances. This behaviour have drastic consequences on the kernels of the integral equations dealing with metallic gratings, which are directly derived from the Green function: the two variable functions relative to the metallic part of the grating tend to delta functions as well. The integration of such functions through a point matching method requires more and more points of discretization around s' = s and, since s can take any value in the interval (0, s d ) the integration and the inversion of the final linear system of equations (bearing in mind that its size is the total number of discretization points) becomes impossible. This remarks explains why the first attempts at implementing the integral equations on computers were not able to give any result for metallic gratings in the visible and infrared regions.

A very efficient way to overcome this difficulty is to apply an approach called local summation [START_REF] Maystre | Thèse d'Etat[END_REF], using another form of the Green function [START_REF] Abramowitz | Handbook of mathematical functions[END_REF]:

( )

0 imd 0 m 1 (r r ') e H k r r ' md x , r (x, z) 4i α + - -- = ∑      G = , (4.101)
with x being the unit vector of the x axis. This form is the direct consequence of the fact that ( )

0 1 H k r r ' 4i + -   is the Green function of the Helmholtz equation: ( ) ( ) ( ) 2 2 0 0 1 1 H k r r ' H k r r ' r r ' 4i 4i k + +     - - -        +  ∇ = δ       . (4.102)
Since the pseudo-periodic Green function (r r ')

-   G in vacuum is defined by: m i d 2 2 m , (r, r ') k (r, r ') e (r r ' md x) α =-∞ +∞ ∇ + = δ -- ∑         G G , (4.103) it follows that (r r ') -   G is a sum of Green functions ( ) 0 1 Ĥ k r r ' md x 4i + --   satisfying: ( ) ( ) ( ) n n n i d i d 2 2 i d 0 0 ˆê r r ' md x k e r r ' md x ˆ e r 1 1 r ' m H k H k 4i d x 4i . + + α α α            ∇ -- + -- = = δ --        (4.104)
Inserting the value of (r r ') -  G given by eq.

∞ ′ α - + -∞ ′ ′ ′ ′ ′ φ - φ ∫ ∫   G = (4.105)
In the same way it can be shown that

d 0 s i (x x) 1 0 dz z z (x x) ik ds (s,s ) (s )ds H (k r r ' ) e (s )ds . 4 r r ' ∞ ′ α - + -∞ ′ ′ -- - ′ ′ ′ ′ ′ ′ ψ - ψ - ∫ ∫     N = (4.106)
When the permittivity becomes very large in modulus, the functions (s,s ) ( )

1 K k n r r ' -   -      
(see [START_REF] Abramowitz | Handbook of mathematical functions[END_REF]) and tend to delta functions when n -→ ∞ . Thus, these functions vary much more rapidly than the unknown function φ , which can be considered as a constant. Thus, remarking that when

x x ′ ≈ [ ] 2 2 2 2 0 dz r r ' x x 1 , dx dz 1 d z z z (x x) (x x) , dx 2 dx exp i (x x) 1,   ′ - - +     ′ ′ ′ --- ≈ - - ′ - ≈ α ≈   (4.107)
yields finally:

( ) d 0 s i (x x) 0 2 0 (s) (s) (s,s ) (s )ds H (k n r r ' ) e ds , 4i dz 2ik 1 dx ∞ ′ α - + - -∞ φ φ ′ ′ ′ ′ φ ≈ - ≈   +     ∫ ∫   G (4.108) and d 2 s 2 3/2 2 0 d z (s) dx (s,s ) (s )ds dz 4ik 1 dx ψ ′ ′ ′ ψ ≈     +           ∫ N . (4.109)
It is not surprising to note that in this approximation, the calculations of the integrals require neither summation of the kernels, nor integration of the singularities, nor matrix multiplication. As the kernels tend toward delta-distributions with amplitudes determining the coefficients in eqs.(4.108) and (4.109), their matrix representations tend to diagonal matrices. Numerical results have shown that this simpler formulation not only successfully applies in the domain where the summation and integration processes defined in the previous sections fail, but also remains valid with a good accuracy (about 10 -3 in relative value) in a large domain of metals and wavelengths. For example, with aluminum, this approach works in the visible, a domain in which the classical method of integration can be used as well (but with a greater computation time), wheras the local summation is necessary for metals in the far-infrared and microwaves domain. It is very important to notice that, using the local summation and assuming that | n | -→ ∞ , it can be shown that the integral equation for metallic gratings described in appendix 4.B tends towards the integral equations obtained for perfectly conducting metals.

Problems of edges and non-analytical profiles

When the grating profile has edges or corners (in 2D case), fundamental difficulties appear. First, the uniqueness of the solution of the electromagnetic field is not ensured. The hypothesis of the integrability of the unknown function u(x) in the integral equation is equivalent to the Meixner condition of integrability [START_REF] Meixner | [END_REF], although it is singular.

The second problem lies in the validity of the boundary condition of electromagnetic field on the grating profile. Indeed, on the edges, the normal and tangential direction on the profile are not defined in a unique manner. Moreover, when establishing these boundary conditions, the demonstration does not work if the surface has edges; However, these warnings are not dramatic. The demonstration of Archimede theorem is also questionable if the object presents edges. Does it exist any doubt about the validity of this theorem?

The third problem in the integral method lies in the process of integration in the vicinity of the edges. The kernels become discontinuous or even meaningless, depending on whether the point of calculation of the integral coincides with the edge point or not, see eq.(4.90) -(4.92). Moreover, the integration can fail due to the eventual singularity of the unknown function on the edge [44].

To overcome the edge problem, there exist several approaches. The most direct one consists in replacing the actual profile z = f(x) with its truncated Fourier representation:

M imKx M m m M f (x) f e =- = ∑  . ( 4

.110)

The new profile has no edges and in most cases of ruled or holographic gratings, it mimics quite well the true profile. Numerical tests have shown a very good convergence of results with respect to the number 2M+1 of Fourier terms and the number P of discretization points, provided the empirical rule P > 4(2M + 1). This method can be used to describe some profiles that are not represented by continuous functions, for example lamellar gratings, provided that M is larger than 10. Unfortanately, there are several important cases of classical gratings and classes of relatively new types of gratings and periodic systems that cannot be treated using this simple approach.

The first problem covers the case of echelle gratings, working in grazing incidence in high and very high orders (see Chapter 1). The fact that the period is some tens or hundreds times larger than the wavelength λ, and that the working facet is quite steep (sometimes going up to 86° groove angle), requires that its geometry is represented in the method with an error smaller than λ/20. Simple but incorrect estimations show that this rule of thumb would be acceptable for numerical treatment: if the period is close to 50 wavelengths, we need about 1000 points of discretization and thus 250 Fourier harmonics, according to the rule P > 4(2M + 1). However, Gibbs phenomenon will significantly modify the form of the working facet, which could be almost vertical. In order to correctly describe the field on this facet, we need to have at least 5 to 10 points per wavelength along the the working facet, rather than along its projection on the x-axis. With a 85° groove angle, the length of this facet is about 11 times its x-projection, in such a way that the number of discretization points must be multiplied by a factor 10. The number of Fourier components in the profile follows the same rule. Another class of problems consists of unconventional geometries, like inverted slope grooves, obtained during groove etching technologies, as seen in Figure .4.5, or rod gratings, shown in Figure .4.6. We have noticed that the problem of vertical segments adds difficulties. For example, the problem of edges cannot be solved any longer by a Fourier expansion of the profile. It exists a possibility to simultaneously solve the two problems by introducing a curvilinear coordinate that follows the grating profile. As far as the integration is made along the profile, this curvilinear integration comes as a natural way of calculating the integrals in the integral equations. Moreover, adaptive meshing can be used to reduce the influence of edges.

In the real life, edges do not exist. On each edge there is at least one atom that has no edges, even though the light has a wavelength much larger than the atom dimensions. Edges of nuclear particles are not discussed even in the most exotic theories. Thus the idea is to replace the edges by arcs, where the partial derivatives can be well defined till the second order, which is sufficient for integral equations, as remarked in section 6.2. An adaptive density of discretization points gives the possibility to significantly increase the density of the points close to the initial edges, and not elsewhere. Let us, for example, consider a rod grating having a straight rectangular cross section. The segments 1 and 3 are parallel to Ox, the segments 2 and 4, to Oz. Let us assume the origin of the curvilinear coordinate at the bottomleft corner (figure 4.7). The coordinate s starts at x = R and z = 0, and follows the interface along its different parts: 1) the first horizontal segment, 0 s a 2R : x (a R) z R R --+ -= .

< < - x(s) R s, z (s) 0.
(4.112)

Here,

1,max 0 s s R : 2 π < - < 1,max 1,max 1,max s R s , s s x a R R sin( ) a R R sin , R s s z R R cos( ) R R cos . R = ϕ + - = -+ ϕ = -+ - = - ϕ = -
C = (a -R, R) (a, R)
and the same for the rest of the profile. The process is straightforward and needs an adapted application for each class of profiles. The advantage is that the derivatives dx ds and dz ds exist and are continuous everywhere on the profile. Thus, the second-order derivatives, which are required for the explicit summation and integration of the kernels, exist at least piesewisely, too. This could easily be checked at the point (a -R, 0), for example. For

- = = = - = = = (4.116)
The comparison of (4.115) and (4.116) points out the existence and continuity of the first derivatives.

It is worth noting that, in contrast with the adaptive spatial resolution used in several other methods (FEM, FDTD, RCW, coordinate transformation methods), here it is more convenient to call it adaptive profile resolution method, as far as it represents a 1D curvilinear coordinate adaptation.

Let us impose the requirement that the arc segments require N arc times larger density points than on the straight segments. The entire length along s of the profile in figure.4.7 (left) is equal to L tot = 2a + 2b -8R + 2πR. The total number of discretization points is related to the length ot the segments, R, and to N arc . If the distance between the points along the straight segments is ∆, it will be equal to arc / N ∆ on the arcs. Thus the total number of points for a single-rod per period is equal to:

arc 2a 2b -8R 2 RN P + + π = ∆ . (4.117)
In practice, unless some automatic scheme of determining the technical parameters of the computation is used, this equation determines ∆ when the total number of integration points is chosen.

Starting from s = 0 in figure 4.8 (right), the abscissa values along s of the points on the interface are given by consecutively adding the values of the point number j to the end values of the previous segment: 1) first horizontal segment, s s - = + ∆ , 4) so on to close the rod, then eventually going to another object inside the same grating period.

An additional improvement can be performed by making a smooth transfer from ∆ to arc N ∆ , i.e., to smoothly go from the density defined on the straight segments and on the arcs. This could be important for large-period systems with respect to the wavelength, if there is a restriction in the total number of points in the profile discretization. If so, we can lay on the fact that the smaller the curvature of the segment (i.e., the larger its length), the smoother the behaviour of the kernels and of the unknown functions. Thus in the middle of a straight segment we will place points that are more distant to each other than close to the extremeties of the segments. Maystre has adapted such approach in his code Grating 2000, by defining a specific distance along the large segments starting from their edges, so that the density of the discretization points increases when approaching the ends of the segments. This approach gave the possibility to model echelle gratings working in very high diffraction orders (600 or more) in the '90s, when no alternative approach was able to provide reliable results. The method was unbeatable for echelles, before Li and Chandezon [47] formulated an improvement of the coordinate transformation method to work for profiles with edges. However, the latter does not apply to rod gratings, or to profiles having the form as given in figure.4.5.

Examples of numerical results

All the results shown in this section are obtained using the code Grating 2000.

Sinusoidal perfectly conducting grating

Table 1 shows the efficiencies in the two non-evanescent orders (-1 st and 0 th ) of a sinusoidal perfectly conducting grating of period 600 nm and height 180 nm (from the bottom to the top of the groove) illuminated under incidence angle 30° and wavelength 600 nm. In this case, the -1 st order is scattered with an angle of scattering (measured anticlockwise, in contrast with incident angle) of -30°, which entails that it propagates just in the direction which is the opposite to that of the incident wave (this is called Littrow mounting by specialists of gratings).

The number of discretisation points is P and the series included in the kernel are summed from -M to +M. The symbol It is worth noting that a precision better than 0.01 is reached as soon as P > 4 and M > 2!. A precision of 10 -3 needs P > 50 and M > 20.

Echelette perfectly conducting grating

The echelette grating is a grating with triangular groove (figure 4.8). The blaze angle b (angle of the large facet with the x axis) is equal to 30° and the apex angle A (angle between the two facets) to 90°. The other parameters are the same as in section 7.1. It must be noticed that the incidence angle and the blaze angle are equal, which entails that the incident wavevector is orthogonal to the large facet. In these conditions, it can be shown that for TM polarization, the grating problem can be solved in closed form: the efficiency in the -1 st order is equal to unity while in the 0 th order it vanishes [48].

The demonstration is straightforward: the sum of the incident wave and of a plane wave with unit amplitude propagating in the opposite direction satisfies all the conditions of the boundary value problem stated in section 3.2. The reader can notice that this sum satisfies the Helmholtz equation. In addition, this sum of two plane waves propagating in opposite directions constitute an interference system that presents white areas and dark lines. The maximum of white lines coincide with the large facets of the grating, and on these lines, the derivative of the field (thus the normal derivative) vanishes. The normal derivative with respect to the small facet vanishes as well since the field is invariant in the normal direction. At the first glance, this property is obvious since the field is "reflected" by the large facets, or in other words, the scattering phenomenon reduces in that case to a simple reflection phenomenon. This reasoning fails since the same phenomenon is not observed for TE polarization: it is dangerous to invoke reflection phenomena on the large facets when the width of these facets has the same order of magnitude as the wavelength of the light! The concentration of the incident energy in a single order is called 'blazing effect' and the theorem of Marechal and Stroke is the origin of the name 'blazed gratings' given sometimes to ruled gratings. Table 4.2 shows that the convergence is significantly slower than for the sinusoidal grating, due to the edges. Nevertheless, a precision of about 0.01 is obtained when P > 25 and M > 10. With the same values, the Marechal and Stroke theorem is satisfied with a precision better than 0.003. The computation time is always less than 1second on a PC computer except for the last line, for which 2 seconds are required.

Lamellar perfectly conducting grating

The profile of a lamellar grating is shown in figure 4.9. The widths of the hole and of the bump are denoted by t and b and the height by h. In this example, b = t = 300 nm, and h = 180 nm. The main conclusion to draw from Table 4.3 is that the convergence for lamellar gratings is even slower than for echelette gratings. This is not surprising if we notice that the number of edges is multiplied by 2. A precision of about 0.02 is obtained when P > 50 and M > 20. The results for P = 6 and M = 3 are aberrant, specially for TM polarization.

Aluminum sinusoidal grating in the near infrared

We consider a sinusoidal aluminum grating with period d = 400 nm, a height h = 100 nm, illuminated with incidence angle 10° and wavelength 300 nm. With these parameters, three plane waves (-1 st , 0 th and +1 st ) are reflected. The optical index of aluminum at 300 nm is equal to 4.2+ i 21.5. We give in Table 4.4 the efficiency in the -1 st order and the sum of the three efficiencies. 4.4 shows a very good convergence of the results, similar to the convergence observed for sinusoidal perfectly conducting gratings, thanks to the local summation of the two variable functions of the kernel derived from the Green function in aluminum.

TE polarization TM polarization P,M

-1 ρ m ρ ∑ -1 ρ m ρ ∑

Buried echelette silver grating in the visible.

The buried silver grating is shown in figure 4.10. A symmetric echelette silver grating with period 900 nm has been covered by a dielectric of index 1.5, the maximum depth e of dielectric being equal to 800 nm. This grating is illuminated in normal incidence .by a plane wave of wavelength 600 nm. The index of siver for this wavelength is equal to 0.006 + i 3.75. Table 4.5 gives the efficiency in the 0 th order and the total of efficiencies of the three reflected orders (-1 st , 0 th and +1 st ). The convergence is slower than in the preceding case and the results for P = 10 and M = 4 are not correct, specially for TM polarization. For TE polarization, a convergence of the results with a precision better than 0.006 is obtained for P = 30 and M = 13. This is not so for TM polarization, in which it is necessary to reach P = 100 and M = 45 to get a precision of the order of 0.003. Here, the method of local summation is not used, but this fact does not explain the slower convergence since the modulus of the optical index is not large. The main reason can be found in the edges of the profile. The slower convergence for TM polarization is rather general for metallic gratings, due to the existence of plasmon resonances on the grating surface.

Dielectric rod grating.

The grating is made of dielectric elliptic rods with optical index 1.4, width w = 600 nm and height h = 400 nm (figure 4.11). The period is equal to 800 nm. It is illuminated with incidence 20° by a plane wave with wavelength 600 nm. Two orders (-1 st and 0 th ) are reflected and transmitted. We give in Table 4.6 the efficiency in the 0 th transmitted order and the sum of efficiencies of the 4 scattered orders, which should be equal to 1 for a perfect energy balance. The precision for P = 30 and M = 13 is equal to 0.01 and for P = 100 and M = 90, it reaches 0.0004.

Flat perfectly conducting rod grating

The grating is similar to that of figure 4.11, but its height is very small (8 nm). In order to obtain a significant reflection, the dielectric has been replaced by a perfectly conducting material The convergence is very slow and it is necessary to reach P = 300 and M = 140 to obtain a precision better than 0.006, the energy balance being satisfied with a precision better than 10 -3 . The reason must be found in the small height of the rods, as explained in section 6.1. The functions included in the kernel of the integral equation become very large in modulus and very small in width for two points located on both sides of the rod for the same abscissa, thus the integration needs a large density of discretization points.

TE polarization TM polarization P,M 0 τ m ρ ∑ + m τ ∑ 0 τ m ρ ∑ + m τ ∑ 100,
F x,z in V , F x,z = F x,z in V , - -     
which satisfies the following conditions:

• it is pseudo-periodic along the x axis:

( ) ( ) ( ) 0 F x+d,z F x,z exp iα d , = (4.118)
• it satisfies a Helmholtz equation:

( ) 2 2 2 F k n F 0 in V , ± ± ± ± ∇ + = (4.119)
• it satisfies a radiation condition for z → ±∞ . The aim of this appendix is to use the second Green's identity and basic theorems on boundary value problems in order to find an integral expression of this function and to analyze the properties of this integral expression. We will deduce the basic keys for writing an integral equation from a boundary value problem.

4.A.2. Calculation of the Green function

The first step of the calculation is to find the pseudo-periodic elementary solutions G , which satisfy the radiation conditions for z → ±∞ and the following equations :

( ) m i d 2 2 2 + - m , (r) k (r)
e (r md x), in V and , n V

α + + =-∞ +∞ + ∇ + = δ - ∑    G G (4.120) ( ) m i d 2 2 2 + - m , (r) k (r) e
(r md x), in V and , n V α - - =-∞ +∞ - ∇ + = δ - ∑    G G (4.121) 0 i d (r d x) (r)e α ± ± + =   G G , (4.122)
with:

m 0 2 mK, K d π α = α + = , (4.123)
and x being the unit vector of the x axis. We must emphasize that, in contrast with eq. (4.119), in which F + and F -do not satisfy the same Helmholtz equation, each Green function satisfies a unique Helmholtz equation in the entire space, with constant On the other hand, the right-hand member of eq. ( 4.120) , called Dirac comb, can also be expanded in series: ( )

m m i d i x m , m 1 ê (r ndx) (z) e d α α =-∞ +∞ δ - = δ ∑ ∑  . ( 4
2 m m m 1 G ''(z) G (z) (z) d ± ± ± + γ = δ , (4.126)
with ( ) ( )

2 2 2 2 m m k n ± ± γ = -α .
For z 0 ≠ , eq.( 4.126) becomes the well-known one-dimension propagation equation in a homogeneous media without sources and have for solutions exponentials. We are searching for plane waves that satisfy the radiation condition for z → ±∞ , thus:

m m m m m A exp i z , if z 0, G (z) B exp i z , if z 0. ± ± ± ± ±    γ >    =     -γ <    (4.127)
Distribution theory proves that the solution m G ± of eq.(4.126) is a continuous function of z and that its derivative has a jump equal to 1/d at z 0 = . These two conditions allow one to find the unknown amplitudes:

( ) m m m m m m A B , 1 i A i B , d ± ± ± ± ± ± = γ --γ = (4.128)
and thus:

( )

m m m 1 G (z) exp i z 2i d ± ± ± = γ γ . ( 4.129)
The final form of ± G becomes:

m m m m 1 1 (r) exp i x i z 2id ± ± ±   = α + γ   γ ∑  G . (4.130)
Once the source is not in the origin, the Green functions is the function (r r ')

± -   G .
Notice that ± G symbolizes the Green functions. ( Similar considerations apply to F (x, z) + , so that, after elementary calculations:

S d H H 0 s S i z z imK x x ' m m i z z imKx m m H,m H,m 1 e ( d x x ', z z ' dF (x ', z ') x x ', z z ' F (x ', z ') ds ' dN dN 1 e dx ' 2d 1 s ') (s ') i 1 i e , 2 - - - - - - ∞ γ - - - =-∞ ∞ γ - - - - = - + - -∞   - - - -   φ + ψ       φ + ψ -       = γ   =     γ   ∫ ∑ ∑ ∫ G G (4. ) ( ) ( ) ( ) ( ) m H m H H
( ) ( ) d 0 0 0 0 s i x x ' i x x ' i x ' i x ' s' 0 F (x, z) (x, z,s ')e (s ')e (x, z,s ')e (s ')e ds ' α - α - α α ± ± ± ± ± =   = ± φ + ψ     ∫ G N (4.136) with: ( ) ( ) ( ) ( ) ( ) m imK x x ' s' i z z' s' m m 1 1 (x, z,s ') (x, z, x s ' , z s ' ) e , 2id ± ∞ - + γ - ± ± ± =-∞ = = γ ∑ G G (4.137) ( ) ( ) ( ) ( ) m s imK x x ' i z z ' m m m (x, z,s ') (x, z, x ' s ' , z ' s ' ) N 1 dx ' dz ' sgn z z ' e , 2d
ds ' ds ' U (x, z) (x, z,s ') (s ') (x, z,s ') (s ') ds '

± ± ± ∞ -+ γ - ± =-∞ ∂ = - = ∂   α = --   γ     ∑ G N (4.
± ± ± ± ± ′=   = ± φ + ψ   ∫ G N . ( 4
.139)

4.A.4. Equation of compatibility

In this section, we establish a crucial property of the integral theory of gratings, which unfortunately is ignored in most of the reference books of Electromagnetics. With this aim, it is necessary to point out a fundamental property of the expression of U (x, z) ± given by eqs. The conclusion to draw from these properties is that the limit of F (x, z) This result is not surprising for the specialist of boundary value problems. Indeed, using the second Green theorem in V -, we introduce in the integral expression of the field the limit values of both the field and of its normal derivative below S. In the domain of boundary value problems, it is well known that if a function must satisfy a Helmholtz equation in V -and a radiation condition for z → -∞ , one cannot impose the limit values of both this function and its normal derivative on S. In fact, we can impose either the limit values of this function or that of its normal derivative on S: in both cases, the solution of the boundary value problem exists and is unique. Unfortunately, it does not exist any tool of applied mathematics which enables one to express the field in an integral form including either the limit values of the field, or that of its normal derivative on S: a prior solution of an integral equation is required, which is much more difficult.

On the other hand, if both the actual values of the field and of its normal derivative are known, the field can directly be expressed in an integral form through the second Green theorem, without any integral equation, and this is why the second Green theorem is considered as a basic tool of the integral mathods of scattering. It must be emphasized that in that case, we know a priori that the the limit values of the field and of its normal derivative are the actual ones, thus that they are compatible, which is not the case for arbitrarily chosen limits.

This fundamental property shows that when the second Green theorem is used with unknown values of the limits of the field and of its normal derivative, we must impose to these limits an equation of compatibility. This compatibility is satisfied if we impose to the

-            or ext s dF (x, z) lim dN - -            .
We can use a first equation:

{ } { } 0 i x + ext -ext lim F (x, z) lim F (x, z) (s)e α - - - - = -ψ    , (4.140)
or equivalently:

{ } { } + ext - ext lim U (x, z) lim U (x, z) (s) - - - - = -ψ    . ( 4

.141)

In order to find a second equation, we can consider eq. ( 4.139) which gives the expression of 0 i x U (x, z) F (x, z)e -α --= . If z is fixed, this expression is a Fourier series in x, which is discontinuous on S. It is well known that the value on S of this Fourier series is the average value of the limits on both sides of S, thus:

{ } { } d s + ext - ext s' 0 lim U (x, z) +lim U (x, z) 2 (s,s ') (s ') (s,s ') (s ') ds ', - - - - - - =   = - φ + ψ   ∫     G N (
(s) (s,s ') (s ') (s,s ') (s ') ds ' 0 2 - - - - - = ψ   φ + ψ + =   ∫    G N . ( 4

.145)

Achieving a similar calculation for V + yields the general equation of compatibility in V ± : sgn z s z ' s ' for s=s'. In fact, a careful analysis of this function around s = s' shows that it is continuous and that its limit when s ' s → can be expressed in closed form [44,[START_REF] Maystre | Thèse d'Etat[END_REF], as stated in section 6.3.

d s s' 0 (s) (s,s ') (s ') (s,s ') (s ') ds ' 0 2 ± ± ± ± ± = ψ   φ + ψ =   ∫     G N , (4.146) ( ) ( ) ( ) ( ) ( ) m imK x s x ' s' i z s z' s' m m 1 1 (s,s ') e 2id ± ∞ - + γ - ± ± =-∞ = γ ∑ G , (4.147) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) m imK x s x ' s' i z s z' s' m m m 1 dx ' dz ' (s,s ') sgn z s z ' s ' e . 2d ds ' ds ' ± ∞ - + γ - ± ± =-∞   α = - -   γ     ∑ N ( 4

4.A.5. Generalized compatibility

In the physical problem, the total field in V + includes the incident field. Here, we give an extension of the compatibility equation to the total field which allows a significant simplification of the use of integral theory.

+ T+ + i + i =F =F +F in V , = =F in V . - -  Γ  Γ  Γ  
F + being the field scattered in V + by a grating illuminated by the incident field i F .

- Γ is the expression of the incident field in V -. Thus, in contrast with the preceding section, the function considered in this section does not satisfy a radiation condition at infinity in V + . According to eq. (4.136), the value of F + is given by: ( )

( ) d 0 0 0 0 s i x x ' i x ' s' 0 i x x ' i x ' F (x, z) (x,
α - α + + + = α - α + +  = Φ -φ    + Ψ -ψ   ∫ G N (4.151)
Thus, the integral expression of the total field + Γ is given by:

( ) ( ) ( ) ( ) d 0 0 0 0 s i x x ' i x ' i i s' 0 i x x ' i x ' i (x, z) F (x, z) (x, z,s ')e (s ') (s ') e (x, z,s ')e (s ') (s ') e ds ', α - α + + + = α - α + +  Γ = + Φ -φ    + Ψ -ψ   ∫ G N (4.152)
or, gathering the terms containing the incident field:

( ) ( ) ( ) ( ) d 0 0 0 0 d 0 0 0 0 s i x x ' i x x ' i x ' i x ' i i i s' 0 s i x x ' i x x ' i x ' i x ' s' 0 (x, z) F (x, z) (x,
α + + = α - α - α α + + + + = Γ =   - φ + ψ     + Φ + Ψ     ∫ ∫ G N G N (4.153)
In order to simplify this equation, we consider the part of the integral containing the incident field, the middle line in the equation above. This part is equal to 0. Indeed, we know that

  φ + ψ =   ∫ G N . ( 4

.154)

Thus finally the expression of + Γ is given by: ( )

( ) d 0 0 0 0 s i x x ' i x ' i s' 0 i x x ' i x '
(x, z) F (x, z) (x, z,s ')e (s ')e (x, z,s ')e (s ')e ds '.

α - α + + + = α - α + +  Γ = + Φ    + Ψ   ∫ G N (4.155)
The result is that the expression of the total field in V + can be obtained from its limit value + Ψ on S and from its normal derivative + Φ , by adding the incident field to the integral expression deduced from the Green theorem.

The generalized compatibility equation is derived from the compatibility equations for

F + in V + and for i F in V -: d s s' 0 (s) (s,s ') (s ') (s,s ') (s ') ds ' 0, 2 + + + + + = ψ   φ + ψ - =   ∫ G N (4.156) d s i i i s' 0 (s) (s,s ') (s ') (s,s ') (s ') ds ' 0. 2 + + = ψ   φ + ψ + =   ∫ G N (4.157)
By adding these two equations, we deduce that:

d s i s' 0 (s) (s,s ') (s ') (s,s ') (s ') ds ' 2 + + + + + = Ψ   Φ + Ψ + ψ =   ∫ G N . (4.158)
The remarkable result is that the compatibility equation given in the preceding section can be extended to the total field: the left-hand side of eq. (4.158) represents the expression of the total field on S and the right-hand, one half of its limit on S. Thus the generalized compatibility condition can be stated in the following way:

The compatibility condition for a field in V + including incident and/or diffracted waves and satisfying the two conditions:

• it is pseudo-periodic along the x axis:

( ) ( ) ( ) 0 F x+d,z F x,z exp iα d , = (4.159)
• it satisfies a Helmholtz equation:

2 2 2 F k (n ) F 0 in V , ± ± ± ± ∇ + = (4.160)
can be stated in the following way: The value on S of the total field, obtained by adding the incident wave to the integral expression deduced from the second Green theorem (but in which the limits are those of the total field) is equal to the half of its limit value on S.

A similar generalized compatibility condition can be obtained for a total field in V - when an incident wave propagating upward in V -(supposed to contain a lossless dielectric) illuminates the grating surface, but this case is not worth in the frame of this chapter.

4.A.6. Normal derivative of a field continuous on S.

The calculation of the normal derivative of F ± on the grating surface in the general case is difficult. However, this aim can be reached at least in one case: when it is possible to define both F + and F -which satisfy three conditions:

• F is continuous across S, or equivalently ( ) ( )

s s + - ψ = ψ ,
• F + and F -satisfy the same Helmholtz equation, with constant 2 2 k (n )

+ , or equivalently n n + -= .

• F satisfies a radiation condition at infinity. Due to the second condition, it seems that this case does not make sense: if n n + -= , one cannot expect any scattering phenomenon. However, the study of this purely mathematical problem is crucial, for example in the study of perfectly conducting gratings.

First, it is worth noting that, thanks to the first condition, the gradient of F can be calculated without any use of distributions, which is not the case if F is discontinuous on S. , the expression of F in the entire space is given by adding the right-hand sides of eqs. (4.161) and (4.162):

α α + + + + + =   = φ + ψ   ∫ G N , (4.161)
( ) d 0 s i x s' 0 F (x, z) F (x, z) F (x, z) (x, z,s ')e (s ') (s ') ds ' α + - + + - = = + = φ -φ ∫ G . ( 4

.163)

Thanks to the continuity of F on S, the expression of its value on the profile does not make problem and is given by: ( )

d 0 s i x s' 0 F (s) (s,s ')e (s ') (s ') ds ' α + + - = = φ -φ ∫ G . ( 4

.164)

In order to obtain the normal derivative of F, let us calculate its gradient:

( ) d 0 s i x (x,z) s' 0 F (x, z) (x, z,s ')e (s ') (s ') ds ' α + + - =   ∇ = ∇ φ -φ   ∫ G , ( 4
α + =-∞ +∞ α + -+ γ - α =-∞ +∞ α     γ     ∇ =           ∑ ∑ G . ( 4

.166)

The components of the normal S N  are given by: S dy ds N dx ds ds ds given by eq. (4.166) are Fourier series in x. Using again the property of discontinuous Fourier series on the discontinuity, it can be derived that: ± N , (s,s ') K is continuous and its limit when s ' s → can be expressed in closed form [44,[START_REF] Maystre | Thèse d'Etat[END_REF], as stated in section 6.3.

  -   =          , ( 4
+ ± α + -+ γ - α ± + =-∞ +∞     =         α   - -     γ         ∑ G (

4.A.7. Limit values of a field with continuous normal derivative on S.

We consider a function F satisfying the following conditions

• F has the same normal derivative on both sides of S, or equivalently ( ) ( )

s s + - φ = φ ,
• F + and F -satisfy the same Helmholtz equation, or equivalently n n + -= .

• F satisfies a radiation condition at infinity. The aim of this section is to calculate the limits of F on both parts of S.

From the secong Green theorem (eq. (4.136)), F can be expressed from the values of ( ) • it identifies to the actual physical solution of the grating problem in + V :

s + φ , ( ) s - φ , ( ) s +
F + + Γ ≡  , (4.191)
• it is continuous across S,

• it satisfies a radiation condition for y → ±∞ .

We denote by It is worth noting that the values of the limits of the field and its normal derivative on S that are the two unknown functions in the classical integral theory (section 2), are now expressed from the single function ( ) s Φ . This is not surprising since the limits on S of the function Γ  given by eq. (4.193) satisfy automatically a relation of compatibility, whatever the function ( ) s ' Φ introduced in the integral may be: it is the field generated by a surface current on S. As a consequence, we have not to include in the theory a relation of compatibility in V + , which was the first integral equation in the classical formalism.

4.B.3. Integral equation

The single integral equation will be obtained by writing the relation of compatibility in V -, considered now to be filled by the actual grating material, i.e. a material of index n -. With this aim, we calculate the limits in V -of the field and its normal derivative on S. inserting in the continuity conditions of the field given by eqs (4.8) and (4.10) the the limit values given by eqs (4.195) with q ± given by eq. (4.12).

The field in V -can be deduced from eqs. s ' q (s,s ') s ' (s ',s '') s '' ds '' 2 q (s) 1 (s,s ') (s ') (s ',s '') s '' ds '' ds ' (s,s ') s ' ds ' 0, 2 2

+ - - = = - + + = =    Φ    + φ + Φ           ψ    + ψ + Φ + + Φ =        ∫ ∫ ∫ ∫ G K N G G (4.202)
which yields, after simplification: i i q q (s) 2 2 2 q q It is interesting to notice that this method has been extended to other problems (including 3D problems of scattering) by many specialists of theoretical physics and applied mathematics [49][50][51][52].

+ + - + - - - - -         + + Φ = - + ψ + φ                      N G G K + N G , ( 4

Figure 4 .

 4 Figure 4.1. Notations

±

  s ± the limit values of U ± on S and being the normal derivative of F ± on S.

  being the value of the periodized incident field, obtained from equation (4.5):

  unknown function Φ . Now, we abandon the field generated by the fictitious surface current density Φ , except the integral expressions of (s) and we come back to the actual physical grating problem. The continuity conditions for the tangential components of the field permit the calculation of (

  denoting the set of non

Figure 4 . 2 .

 42 Figure 4.2. Single layer inside a stack of a multilayer grating

Figure. 4 . 3 .

 43 Figure.4.3. Layer with two profiles that are separable by a plane layer

Figure 4 . 4 .

 44 Figure 4.4. Parameters of the incident wave in conical mount.The angle ϕ denotes the angle between the incident wavevector i k  and its projection i xz k 

  .70) with P being the number of matching points. It is worth noting that the value of jp W may differ from the value of ( ) j p W s ,s obtained through the rectangular rule of integration

  exact integral is equal to the 0 th term in the Fourier expansion of a(s):

  us at first evaluate the asymptotic values of m γ and m α for large values of m:

  interval (0, d), because d s s′ -never vanishes when s,s (0, d) ′∈

  finally gathering the infinite sum of integrals on one period into a single integral from -∞ to +∞ yields:

  replacing k by kn -in (s,s ) ′ G and (s,s ) ′ N. Since the value of n -is close to an imaginary number in the visible and infrared regions for usual metals (for example, n -= 1.3 +i 7.11 for aluminum at 650 nm), the Hankel functions 0 H (k r r ' )

Figure 4 . 5 .

 45 Figure 4.5. Schematic representation of an etched grating profile with non-analytical function description.

Figure 4 . 6 .

 46 Figure 4.6. Grating made of inclined rectangular cross-section rods.

Figure 4 . 7 .

 47 Figure 4.7. Left: rounding of corners of a rectangular cross-section rod with side lengths a and b, together with x, y, and s coordinate lines. Right: schematic representation of different straight and arc segments to obtain the links between the Cartesian and the curvilinear coordinates.

  circular rounding of the bottom-right corner, defined by the equation:

Figure 4

 4 Figure 4.8: A ruled grating

Figure 4

 4 Figure 4.9. A lamellar grating

Figure 4 .

 4 Figure 4.10. A symmetric silver grating coated by dielectric

Figure 4

 4 Figure 4.11. A dielectric rod grating

  .125) Introducing equations (4.124) and (4.125) in equation (4.120), multiplying by 0 i x e -α then identifying the coefficients of the Fourier series yield:

  , we can take advantage of the fact that the segment is parallel to the x axis and lies outside the groove region, which entails that H,m φ and H,m ψ are related through the plane wave expansion valid in the homogeneous region (see chapter 2): us to cancel the integral along H S and thus the values of F (x, z)can be determined by an integral along a single groove of S: curvilinear abscissa of the point of S of abscissa d, the origin of curvilinear abscissa being the origin of the Cartesian coordinates system.

( 4 . 1 .

 41 13) and (4.139). Let us suppose that we introduce in eq. (4.136) arbitrary periodic functions (s ') , we cannot expect to obtain the actual value of F (x, z) important, if the point P (x, z) of V -tends towards a point P S of curvilinear abscissa s located on the profile (figure 4.12), the limit of F (x, z) this surprising property, it is necessary to give two results which can be demonstrated using the theory of distributions. Let us give these two fundamental results without demonstration. The integral expression of F (x, z) - inside V -satisfies the Helmholtz equation in Vlike the actual field. What happens to the same integral expression of F (x, z) - , but now calculated in region V + ? Denoting by ext F (x, z) - the function equal to F (x, z) - in V -and tothis integral extension in V + , it is easy to verify that ext F (x, z) normal derivative across S (difference between the values in V + and in V -) are respectively equal to



  vanishes throughout V + . Indeed, if this property is satisfied, the jumps are nothing but in the integral expression are equal to the actual physical values (s ')

  zero, these two conditions being equivalent. Indeed, if{ } + ext lim F (x, z) 0 -= , the expression of ext F (x, z) - in V + satisfies a Helmholtz equation, a radiation condition at infinity and its limit on S vanishes. The obvious solution of this boundary value problem is ext F (x, z) 0 -=  in V + , and we have seen that the solution of this boundary value problem is unique, thus it is the solution. The consequence is that

  of coordinates x,z becomes a point of S of curvilinear abscissa s: eqs. (4.141) and (4.142), we deduce the limits of U -  on both sides of S, and we derive the equation of compatibility in V -

Figure 4 .

 4 Figure 4.13. generalized: generalized compatibility

  is to be noticed that, when z is constant, the components of

  ordinate of the top of the grating profile. It must be noticed that a finite number of orders m, called propagating orders, are non-evanescent and propagate at infinity. They correspond to real values of m + γ (for reflected orders) or m γ (for transmitted orders, if the optical index n -is real only).Equation (4.184) can easily be generalized to the case in which we know the limit value of the total field on S (including incident waves) and its normal derivative. It suffices to analyze the behaviour at infinity of eq. (4.155) instead of eq. (4.136). The result is that it suffices to replace the values (s) ρ in the reflected orders propagating above the grating can be obtained by using the Poynting theorem on segments of one period parallel to the x axis: is made of a lossless dielectric, the transmitted efficiencies m

  of Γ  and of its normal derivative S d dN Γ  on S, bearing in mind that by definition, ψ ' ψ + + ≡ and ' + + φ ≡ φ .

  (4.199) and (4.200) using eq. (4.136): z,s ') (s ') (x, z,s ') (s ') ds ' of compatibility is given by eq. (4.146):

  -, in such a way that the expression of F + at infinity becomes a sum of plane waves: diffraction efficiencies of the reflected waves are then deduced by: can be derived from eq. (4.201) that the transmitted field can be represented by a sum of plane waves below the profile: the transmitted plane waves are derived from eq. (4.185) after calculating the functions - φ and - ψ from eqs. (4.199) and (4.200):

  The amplitudes m r and m t of the plane waves reflected and transmitted by the grating can be deduced from the solution of the integral equation using eqs. (4.184) and (4.185) of appendix 4.A:

	m r	1 2d	d s 0 s	e	( ) imKx s i z s ( ) m	m i (s)	dx(s) ds	m m	dz(s) ds	(s) ds,

23)

Eqs

(4.19

) and (4.23) constitute a system of two integral equations with two unknown functions, which can be solved on a computer.

.36) 4.3.2. Perfectly conducting gratings for TM polarization

	Once again, the generalized compatibility condition is used (eq. (4.158):
	d s' 0 s	(s,s ') (s ')	(s,s ') (s ') ds '	i	( ) s	(s) 2	.

Table 4 .

 4 1. The sinusoidal perfectly conducting grating.

Table 4 .

 4 2. The echelette perfectly conducting grating and the Marechal and Stroke theorem.

	6,3	0.6531	0.4451	1.0982	1.4452	0.0012	1.4464
	25,10	0.5838	0.4123	0.9961	0.9976	0.0001	0.9977
	50,20	0.5838	0.4123	0.9961	0.9976	0.0001	0.9977
	110,50	0.5929	0.4055	0.9984	0.9984	0.0000	0.9984
	250,100	0.5932	0.4035	0.9967	0.9989	0.0000	0.9989

Table 4 .

 4

	6,3	0.7770	0.6999	1.4769	23.66	5.09	28.75
	25,10	0.3490	0.6394	0.9884	0.8293	0.1724	1.0016
	50,20	0.3347	0.6569	0.9915	0.8191	0.1705	0.9896
	150,70	0.3279	0.6659	0.9938	0.8201	0.1718	0.9919
	250,100	0.3288	0.6733	1.0021	0.8214	0.1725	0.9939

3. The echelette perfectly conducting grating and the Marechal and Stroke theorem

Table 4 .

 4

	6,3	0.5205	0.9582	0.4367	0.9618
	10,4	0.5201	0.9649	0.4325	0.9521
	30,13	0.5203	0.9655	0.4321	0.9518
	100,45	0.5204	0.9655	0.4320	0.9518

4. The aluminum sinusoidal grating

Table

Table 4 . 5

 45

	10,4	0.4892D	0.9830	0.6316	1.0622
	30,13	0.5164	0.9651	0.5665	0.9339
	100,45	0.5153	0.9600	0.6062	0.9178
	200,90	0.5153	0.9593	0.6168	0.9153

. The buried silver grating

• it has the same pseudo-periodicity as the actual solution of the grating problem: 0 i d

Appendix 4.A. Mathematical bases of the integral theory

4.A.1. Presentation of the mathematical problem

We consider (figure 4.1) a function ( ) ( ) ( )

4.A.3. Integral expression

Now, we apply the second Green theorem in order to find the expression of the function F ± . First, we consider the expression of F -in V -:

( ) ( ) ( ) ) is opposite to that of ( )

), which entails that given by eq. ( 4.174) vanishes in V -. The same property holds for the expression of F (x, z) given by eq. ( 4.175), which vanishes in V + , thus we can write that, if the compatibility equations are satisfied:

)

In order to express the limits

on both parts of S,. we can use a first equation:

To find a second equation, we can consider eqs. (4.176) and (4.177) which gives the expression of U(x, z) . If z is fixed, the expression of U(x, z) is a Fourier series in x, which is discontinuous on S. The value on S of this Fourier series is the average value of the limits on both sides of S, thus:

.181)

4.A.8. Calculation of the amplitudes of the plane wave expansions at infinity.

For z → ±∞ , the expression of F ± given by eq. ( 4.136) can be simplified since ( )

( ) z z ' z z ' -= ± -, in such a way that the expression of F at infinity becomes a sum of plane waves:

4.B.1. Definition of the unknown function

The single unknown function Φ is defined from a function

satisfying the following conditions:

• it satisfies the same Helmholtz equation in the entire space, except may be on the profile S,:

( )

The question which arises is to know if - Γ  is well defined. The continuity of Γ  across S imposes the value of - Γ  on S: ψ ' ψ ' ψ -+ + ≡ ≡ . In addition, - Γ  satisfies a Helmholtz equation and a radiation condition at infinity. It is worth noting that, in contrast with the function Γ introduced in Appendix 4.A, here the function Γ  does not include the incident field and thus has no physical meaning below the profile. As mentioned in Appendix 4.A, the solution of this boundary value problem (since we impose that - + Γ = Γ   ) exists and is unique, thus Γ  is correctly defined in the entire space. The unknown function Φ of the integral equation is defined as the jump of the normal derivative of Γ  across S, more precisely: ' '

As regards the physical interpretation of Φ , it can be shown easily that Φ 0 i x e α is the surface current density which, placed on S, generates in + V the actual scattered field. In other words, we have considered a fictitious structure consisting of an infinitely thin, perfectly conducting metallic sheet supporting a surface current density ( ) ĵ s y placed on the grating surface, this surface separating two media having identical optical properties (refractive index n + ). The unknown Φ is proportional to j.

4.B.2. Expression of the scattered field, its limit on S and its normal derivative from Φ .

It must be noticed that the function Γ  satisfies all the conditions of the function F of section 4.A.5. Since Γ  is continuous on S, the calculation of Γ  from Φ can be achieved using eqs.