Overview of Computational Methods for Photonic Crystals

A photonic crystal (PC) is a periodic structure whose refraction index of the material is periodically modulated on the wavelength scale to affect the electromagnetic wave propagation by creating photonic band gaps. In 1887, Lord Rayleigh is the first to show a band gap in one-dimensional periodic structures i.e. a Bragg mirror. In 1987, Eli Yablonovitch and Sajeev John have extended the band gap concept to the two and threedimensional structures and for the first time, they use the term "photonic crystal" (Yablonovitch, 1987; John, 1987).


Introduction
A photonic crystal (PC) is a periodic structure whose refraction index of the material is periodically modulated on the wavelength scale to affect the electromagnetic wave propagation by creating photonic band gaps.In 1887, Lord Rayleigh is the first to show a band gap in one-dimensional periodic structures i.e. a Bragg mirror.In 1987, Eli Yablonovitch and Sajeev John have extended the band gap concept to the two and threedimensional structures and for the first time, they use the term "photonic crystal" (Yablonovitch, 1987;John, 1987).
Progress in computational methods for the photonic crystals is understood through an historical review (Oyhenart, 2005).At the beginning of research in the photonic crystals, the purpose was to find a structure with complete band gap by improving the computational methods.In 1988, John shows theoretically by the scalar method of Korringa-Kohn-Rostoker (KKR) that the face centered cubic lattice (FCC) has a complete band gap between the second and the third band.One year later, Yablonovitch builds this structure and finds a band gap experimentally but the W-point raises a problem.In 1990, Satpathy et al. and Leung et al. confirm the complete band gap by the scalar plane wave method (PWM).A few months later, these two teams improve their methods to obtain vectorial PWM on D and E fields.They find that FCC structure does not have complete band gap because W-point and U-point are degenerate.With these results, the editor of the journal "Nature" writes "Photonic Crystals bite the dust" (Maddox, 1990).Only two weeks later, Ho et al. created the vectorial PWM on H and they do not find the complete band gap in FCC structure but they show a complete band gap in the diamond lattice.In 1992, Sözuer et al. improve convergence of the PWM and they obtain a complete band gap for FCC lattice between 8th and 9th band.This structure that has caused many discussions has a complete band gap but not where it was expected.(1)


In these laws, the same physical behavior is observed if we change simultaneously the wavelength and the structure dimensions in the same proportions.Therefore, it is convenient to introduce a normalized wavelength 0 /a and a normalized frequency af/c= a/ 0 , with a the lattice constant of the photonic crystal (Joannopoulos et al., 2008).
Some methods do not solve the Maxwell equations directly but they use the Helmholtz equations, for example the E-wave equation or the H-wave equation: (2) Equations number to be solved depends on dimension of the photonic crystal.In the twodimensional case, the problem is simplified.It is assumed that the materials are uniform along z-axis.It follows that the fields are uniform and the partial derivatives with respect to the variable z vanish.The previous equation is simplified and split-up into TE-polarization and TM-polarization; we have a scalar equation for each polarization: The computation time decreases exponentially compared to the three-dimensional case.In the one-dimensional case, the problem is even more simplified.If the materials are uniform along x-axis and y-axis, we solve analytically one second-order equation:

Symmetries of photonic crystals
Like all sets of differential equations, Maxwell's equations cannot be uniquely solved without a suitable set of boundary conditions.Photonic crystals have symmetries which define boundary conditions.Symmetries of the structure are not a sufficient condition to reduce the computational domain; the electromagnetic field must be also symmetrical.On figure 1, we study symmetries on a lattice of cylinders for different polarization of the incident wave.In the first case, we apply a TM incident wave on the lattice with a horizontal symmetry and the tangential magnetic field vanishes on the axis of symmetry.We put perfect magnetic conductor (PMC) as boundary condition on the symmetry axis and we study only half top or bottom of the problem.In the second case, we apply a TE incident wave.Similarly, we reduce the problem with a perfect electric conductor condition (PEC).

Fig. 1. Both kinds of lateral symmetries
On figure 1, the geometry has also a vertical symmetry.However, the electromagnetic field is not symmetrical because the incident wave comes only from the left side.A solution is to divide the incident wave into an even mode and an odd mode (figure 2).For the even mode, we have a vertical symmetry plane of the electric field.Only the left or right part is solved with a perfect magnetic conductor condition on the axis of symmetry.The antisymmetric problem is reduced in the same way with a perfect electric conductor condition.The sum of the symmetric and antisymmetric problem provides the solution of the total problem.The reflection of total problem is the sum r even + r odd and the transmission is the subtraction r evenr odd .Two half-problems are solved more quickly than the whole problem because the computation time increases exponentially with the size of the problem.

Periodicities of photonic crystals
Photonic crystals, like the familiar crystals of atoms have discrete translational symmetry T. The dielectric function is periodic therefore the electric and magnetic fields can be written as the product of a plane wave envelope and a periodic function: The result above is commonly known as Bloch's theorem in solid state physics (Kittel, 2005).
If we apply these conditions to an infinite lattice of cylinders (figure 3a), the calculation of band structure is reduced to the study of a single cylinder.

Fig. 3. Periodic boundary on photonic crystals
Figure 3b is a lattice of cylinders, finite along y-axis and infinite along x-axis for studying the transmission of the structure.If we know the transfer matrix of one layer, it is easy to find the transmission of the total structure thanks to the transfer matrices method (TMM).It reduces the computational domain to one layer.It will be detailed in the next sections.
By mixing TMM along y-axis and Bloch's conditions along x-axis, the structure of figure 3 can be reduced to a single cylinder.Figure 4 summarizes these techniques in the 3D case.

Finite element method (FEM)
The finite element method is very popular in mechanics and civil engineering.It was originally developed in the 40s to solve problems of mechanical structures.A few years later, it has been applied to electromagnetism.Since then, the finite element method extends to all branches of physics and engineering where there exists a partial differential equation (PDE) with boundary conditions.It can be formulated from the variational method or the weighted residual method.We will present in outline the second method which is simpler.
Let us consider a partial differential equation with  the differential operator of order n applied to a function φ and a source function f: The first step is to expand the function φ on a set of functions: 1  are the chosen expansion functions and j c are constant coefficients to be determined.The best solution of the equation 6 is obtained when the residue rf    is weakest on all points of the domain .The weighted residual method requires this condition: If we use the Galerkin method, the weight functions w i are the functions of previous interpolations and we write the problem in matrix form: So now, we have to solve a linear system with c j unknowns where most of the entries of the matrix L are zero.Such matrices are known as sparse matrices, there are efficient solvers for such problems.
The basic idea of the finite element method is to divide the computation domain into small subdomains, which are called finite elements, and then use simple functions, such as linear and quadratic functions, to approximate the unknown solution over each element.For plane geometries, the domain is divided into finite triangular sub-domains.For three-dimensional problems, the sub-domains are tetrahedra.These two-and three-dimensional finite elements are widely used because it is a variable mesh and adapts to curved structures.

Transmission calculation of a photonic crystal
The finite element method usually solves the E-wave equation for PC (equation 5): For most of the finite elements computer programs, the frequency is fixed and the electric field is the unknown (Massaro et al., 2008).We use a commercial software, Ansys HFSS (High Frequency Structure Simulator).This three-dimensional computer program builds an adaptive tetrahedral mesh to model microwave devices, for example microstrips and antennas.These devices have a characteristic length lower than the wavelength of study.It is more difficult to study PC because it has a periodicity close to the wavelength and the matrix of calculation is large.To simplify calculations, we will study a periodic PC according to two directions of space.The directions where PC is infinitely periodic require to study only one period according to this direction.Several methods exist to model these conditions.
In the first method, the source is a plane wave and we use periodic boundary conditions (PBC) on lateral faces and absorbing boundary conditions (ABC) on bases.We study only one lateral period (figure 5).If the source has an oblique incidence, the phase shift is easily taken into account by PBC because the E-field is a complex vector in FEM.The second method requires a very different source, a wave port.This source is a semiinfinite waveguide whose cross section is drawn on bases of structure.Propagative modes of this fictitious guide will be the source of the structure.For a plane wave source, we apply to lateral faces and on the fictitious guide the conditions PEC and PMC (figure 6).The source is transverse and has a normal incidence.It is not possible to change the angle of incidence without change the boundary conditions of the lateral faces.The PC studied on figure 5 can be reduced to quarter-spheres thanks to the symmetry of the fields at normal incidence.
Fig. 6.Calculation of E-field in the PC with PEC, PMC and wave ports

Band structure calculation
The band structure shows the states which propagate in a PC.These states are differentiated by their frequency and their Bloch wave vector.FEM sets the wave vector and solves the wave equation to find the frequencies.There is no source, only boundary conditions to set because it is an eigenvalue problem.In the case of a cavity resonator, the boundary conditions of the domain are PEC.Whereas, for PC, we choose the Bloch conditions on all faces of the unit cell to set the wave vector.The phase shift of the Bloch conditions is set easily because the fields are complex vectors.The FEM calculates the band structure of dielectric with or without losses, metallic, and metallodielectric PC.Any material can be used by this method, it is the main advantage.
In this chapter, we choose a photonic crystal to study, a cubic lattice with several layers of dielectric spheres which have a permittivity equal to 5.1 and radius equal to 0.4*a (a: lattice Source: incident wave constant).For the band structure, the number of layer is infinite, and for transmission calculation, we study four and eight layers.The structure is deliberately simple to be studied by all the methods.The band structure, the transmission coefficient in normal incidence and some modes of dielectric PC are plotted on figure 7.For normal incidence, the first two band gaps are found in the band structure and transmission curves.
The first method with incident waves takes about 7 hours and 960 Mo of working memory on a personal computer to calculate the transmission of 8 layers i.e. 8 spheres.Whereas if we use the second method with the wave ports, the memory is reduced to 360 Mo and the CPU time is reduced to 8 min.To obtain this optimization, we reduced the geometry to four quarter-spheres thanks to symmetries and we use the Padé interpolation on frequencies.It is necessary to make optimizations if you want to use FEM.

Finite-difference time-domain method (FDTD)
Finite difference method is a numerical method for approximating differential equations using finite difference equations to approximate derivatives.This is very simple to implement but it has some limitations in the mesh geometries.In 1966, Yee proposed a finite difference scheme applied to electromagnetism.The FDTD was born (Taflove & Hagness, 2005).The success of this method is due to the scheme based on the Taylor series of second-order.For example, this method is used to model the effects of cellphones on the human body, antennas and printed circuits.The FDTD uses the two structural equations of Maxwell in a conducting, isotropic and homogeneous medium: www.intechopen.comPhotonic Crystals -Introduction, Applications and Theory

274
The FDTD uses an approximation of derivatives by centered finite differences.
Let us apply this approximation to one-dimensional case in order to understand the FDTD principle.The E and H-field are stepped in time and space.By replacing curls and derivatives, we obtain one-dimensional Yee scheme: The E and H-field are staggered and updated step by step in time.E-field updates are conducted midway during each time-step between successive H-field updates, and conversely.This explicit time-stepping scheme avoids the need to solve simultaneous equations, and furthermore it is order N i.e. proportional to the size of the system to model.This scheme can be generalized for two-dimensional and three-dimensional problems.The Yee scheme is stable if the wave propagates from one cell to another with a speed less than the light (the Current-Friedrichs-Lewy condition).

Calculation of transmission coefficient
FDTD is used to study PC.Calculation domain of finite PC is surrounded by absorbing boundary conditions (ABC).The periodic infinite PC uses periodic boundary conditions (PBC) on the lateral faces (Figure 8).The source is Gaussian with a spectrum which extends on the frequency range to study.The fields are calculated on the time domain and we use the Fourier transform to convert them on frequency domain.
Fig. 8. Source, probe and boundary conditions for infinite PC For the non-normal incidence, a simple periodic boundary condition cannot be applied but many methods exist to solve this problem.One method is to calculate the computing problem twice, with the sources cos(t) and sin(t).The addition of two calculations is the solution of the source exp(-it).Thanks to this source, we can apply the Bloch conditions on complex vector fields.This method is simplified for the study of band structure.

Calculation of band structure
The band structure is calculated from the eigenmodes but FDTD is not an eigenvalue problem.We use the unit cell of PC with boundary conditions of Bloch and several Gaussian functions for the source.We only need one calculation to apply the Bloch conditions.After 100000 time-steps of calculation, propagative waves are amplified and evanescent waves vanished.If the number of time steps is too small, the transmission peaks are widened, therefore imprecise.To reduce the number of time-steps without affecting the accuracy, we can use the Padé approximation or ADI-FDTD formulation (Taflove & Hagness, 2005).On Figure 10, we plot the electric field amplitude of the unit cell of PC.This curve is similar to the diffraction pattern in solid state physics.It is plotted for all points of the first Brillouin zone to view the band structure.If we compare the curve of figure 10 with the band structure calculated by the FEM, the mode 0.57 is much attenuated and the following mode does not exist on figure 10.In fact, the source was not correctly selected to excite these modes, and so we must ensure that the source excites all available modes.

Finite-difference frequency-domain method (FDFD)
We will study a finite difference method combined with transfer matrix method and the Fourier transform.The finite difference method substitutes the derivatives in PDE to obtain finite difference schemes.In electromagnetism, FDFD uses the structural equations of Maxwell on space (k, ) and applies approximations on the wave vector: a,b,c are lattice constants along x, y and z-axis.We obtain discrete equations which are applied to one layer of the structure on a cubic lattice.In 1992, Pendry and MacKinnon (Pendry & MacKinnon, 1992) used this method with transfer matrix method (TMM) which extended the solution of one layer to the total structure.To apply TMM, the discrete equations must be written on real space.We obtain a set of equations where z-component of the E and H-fields can be removed.The six equations are reduced to four equations which are written in a matrix form: (15)

Calculation of band structure
The band structure is calculated from unit cell of structure by setting frequency and calculating wave vectors propagating in the unit cell.Bloch's theorem applies to vector F: If the above equations are joined together, we have an eigenvalue problem: We set the wave number k  and the wavelength k 0 to solve the eigenvalue problem.The eigenvalues k z having an imaginary part are eliminated because they are not propagating waves.The remaining eigenvalues gives us the band structure k z (k 0 ).

Calculation of transmission coefficient
The reflection and transmission calculation is performed using the transfer matrices.It is more interesting to calculate the elements of transmission matrix than the transmitted field in some points.The incident, reflected and transmitted waves are expanded on a planewave basis sets thanks to the transfer matrices.The right and left eigenvectors are different because ˆ0 T operator is not self-adjoint.They are expanded on the reciprocal space k x and k y and translate on the direct space with the Fourier transform because the transfer matrix is known in direct space.k z is found easily because we are in a vacuum layer.The transfer matrix T is not compatible with the new plane-wave basis set, we convert this matrix: To calculate the transmission coefficient, this matrix is arranged to group the eigenvectors with the same propagation direction.This matrix is divided into four blocks: TMM calculates the transmission coefficient of several layers from one layer (see 1D-MST).A great number of layers can be calculated without increasing the computing power.
Figure 11 plot the band structure and the transmission coefficient calculated by FDFD, for the previous dielectric PC.If we compares with figure 7, the first bands are correct but inaccuracies on the higher band after 0.6-frequency are due to the weak mesh (7x7x7 cells).
If we increase the number of cells, accuracy increases.Calculation is more difficult for high frequencies because the E-and H-field are functions that oscillate more.The transmission calculation of 8 layers is calculated in 22 seconds with 3.3 Mo of working memory for 7x7x7 cells.The computing time is low because only one sphere is actually calculated thanks TMM.

Finite Integration Technique (FIT)
In 1977, Weiland proposes a spatial discretization scheme to solve the integral equations of Maxwell (Weiland, 1977).This scheme called Finite Integration Technique (FIT) can be applied to many electromagnetism problems, in time and frequency domain, from static up to high frequency.The basic idea of this approach is to apply the Maxwell equations in integral form to a set of staggered grids like FDTD (figure 12).The spatial discretization process is applied to the integral form of the Faraday's law: x y z td td t ei j ke i j ke i j ke i j kb i j k t This last equation is an exact form without approximation.This process is applied to the other Maxwell equations in integral form.CST Microwave Studio is the software based on the FIT.Unlike FDTD, we use the local integral form and we can apply the technique of Perfect Boundary Approximation (PBA) which decreases the meshes on the boundaries.

The tight binding method (TB)
The tight binding method (Lidorikis et al., 1998) is less used than the PWM although this method is fast for the calculation of defect states in a PC.By analogy with the TB model in solid-state physics, H-field is expanded on localized wave functions.They are calculated from the Wannier functions a n (R,r) for a PC without defect:   The Wannier functions are a complete set of orthogonal functions which is defined for each band and each unit cell.The Wannier functions only depend on the quantity (r-R), r is the space position and R is any lattice vector.The reverse relation is written as follows:   H n (k,r) functions are the eigenvectors of H-wave equation.PWM solve this equation to get the H n (k,r) functions for the PC without defect.To study a defect in PC, the H-field is expanded on the previous Wannier functions: Green's functions method solve the H-wave equation for the photonic crystal with defect.

Plane Wave Method (PWM)
Plane wave method is used to study the band structure of PC.It comes from the solid state physics where the electronic wave functions are scalar whereas the electromagnetic fields are vectors.A scalar approximation of fields is not enough to describe band structure.This method is modified to take the vectorial nature of the fields into account.Three computational methods of vectorial PWM were created quasi-simultaneously: the E-field method (Leung & Liu, 1990), the D-field method (Zhang & Satpathy, 1990) and the H-field method (Ho et al., 1990).The D and H-field are continuous in PC unlike the E-field.Moreover, only the differential operators of the wave equation on E-field and H-field are self-adjoint.We present the H-field method because of these properties.PWM expands the field and permittivity on plane-wave basis set.The permittivity is a periodic function with translational symmetry T, so Bloch's theorem can be applied to the H-field: The H-field is written: As the H-field is transverse, each plane wave is perpendicular to the propagation vector K+G.The transverse plane of the propagation vector is described by the unit vectors 1 ˆ KG e and 2 ˆ KG e . The set of vectors   eeK G represents an orthonormal basis.We only need to storage two vectors instead of three, consequently data storage is reduced.The Fourier series expansion is replaced in the H-wave equation: Some algebraic calculations simplify this equation and we get for every vector G, the central equation of the photonic crystals: This equation is an eigenvalue problem which is solved with classical methods.The Bloch vector K is set and we try to find the eigenvalues k 0 .The calculation convergence depends on the N number of reciprocal lattice vectors G.A minimum number of vectors G is necessary to describe correctly the permittivity of the PC.The convergence of the problem is rather slow.As the H-field is transverse, the number of equations is decreased from 3N to 2N.The PWM is difficult to apply to the materials whose permittivity depends on the frequency like metals.On figure 13, we plot the band structure of the dielectric PC.It is calculated in 12 seconds with 10 Mo of working memory for 343 plane-waves.Calculation is fast because the structure is simple.If we compares with figure 7, we obtain the same result.
The 2N equations split into two parts, the TE-polarization and the TM-polarization.On the TM-polarization, the H-field vanishes on the xy-plane and the central equation is written: On the TE-polarization, the H-field is parallel to x-axis and the equation is written: In the one-dimensional case, the equations for TE and TM-polarization become similar.We solve N equations.

One-dimensional multiple-scattering theory (1D-MST)
The multiple-scattering studies the interaction between objects using the analytical solution for each object taken individually.In the next section, we will apply the MST to the cylinders and the spheres.Before, we will establish the analytical solution of onedimensional PC i.e. multilayer structures.Transfer-matrix method is known since many years (Born & Wolf, 1999).It is essential for the study of PC.TMM reduces the computational domain.In this section, transfer matrix will be applied to one-dimensional PC.For oblique incidences, we solve the E-wave equation for TE-polarization and the H-wave equation for TM-polarization.Two polarizations are separated and we use the same steps of calculation for two polarizations.
We will study only the E-wave equation: Let us suppose that the layers of 1D-PC are stacked up along e z .The PC is uniform according to e ┴ and e ║ (figure 15).Because of these symmetries, the E-field and E-wave equation are simplified: We calculate the transfer matrix between layer 1 and 2. The phase shift between the beginning and the end of layer 1 can be written in a transfer matrix form:  The product of previous matrices provides the transfer matrix from layer 1 to layer 2: The transmission of the layer is calculated from inversion of transfer matrix.The inverse will be calculated directly from phase shift matrices and interface matrix to avoid inaccuracies:  Several layers are calculated similarly.

Kronig-Penney model
The Kronig-Penney model evaluates the electronic levels of a crystal structure in a onedimensional periodic potential (Kittel, 2005).This model has been modified to be used in PC and takes into account the oblique incidences (Mishra & Satpathy, 2003).We use TMM and the field expansion on backward and forward wave.The structure is periodic according to zaxis (figure 16).K is the wave vector of Bloch.In the unit cell, the permittivity is not uniform according to zaxis.It is splitted into sub-cells with a constant permittivity.We calculate the transfer matrix of the unit cell from the expressions of the previous section: We expand equation 44 on forward and backward waves: The Bloch factor e iKa and the complex conjugate value e -iKa are the two eigenvalues of the transfer matrix because the determinant of this matrix is equal to one.The trace of the transfer matrix is equal to the sum of the eigenvalues: We get an equation similar to electronic case.The transmission coefficient is different in the TE and TM-polarization.To plot band structure f( k ║ ,K, k 0 ) = 0, we set the wave numbers k 0 and k ║ and calculate the Bloch number K.

Two-dimensional multiple-scattering theory (2D-MST)
The multiple-scattering is an analytical theory which calculates the scattering of N objects from the scattering of each object independently (Felbacq et al., 1994).In the twodimensional case, objects are cylinders and the theory uses the scalar wave equation: Using Green's theorem, the scattered field outside the cylinder can be written as follows: The function H is a Hankel function.The surface integral can be restricted to the Cj cylinders.The scattered field is written as a sum: The cylinders interact to give the scattered field by the structure.To better understand the process of multiple-scattering between the cylinders, we will study a simple example of cylinders with circular section aligned along x-axis and excited by an incident field propagating along x-axis (figure 17).

Fig. 17. Incident field on n cylinders aligned
The incident field is expanded into Bessel functions: We define the total incident field around the cylinder i which takes into account other cylinders: In the case of the cylinders with circular section, the integral 54 is calculated analytically: S jl is the scattering coefficients (Bohren & Huffman, 1998).The field  j s E r scattered by a cylinder j is shifting on another cylinder i to become an incident field G jil .bjl on this cylinder with G jil the translation coefficients (Felbacq et al., 1994).If we shift all the scattered fields on the cylinder i, and if we add the initial incident field to it, we get the total incident field on the cylinder i: G jil and S jl are the translation and scattering matrix coefficients.To get the total incident field and the field scattered by the structure, it is necessary to calculate G jil and S jl coefficients and to solve the multiple-scattering equation.This method is suitable for the study of defects in PC because it does not impose any condition on the position and the material of cylinders.
On figure 18, we study a 2D-PC with 80 cylinders doped by a microcavity.The three-dimensional multiple-scattering theory use same principles as the twodimensional case but three-dimensional case is more complex because the Helmholtz equation is vectorial.The 3D-MST is separated into two parts, the calculation for one sphere and the generalization to N spheres (Oyhenart & Vignéras, 2007).The first part calculates the scattered wave by one sphere with Mie theory.The base of spherical harmonics is written as follows:    The elements of S-matrix are Mie's coefficients (Bohren & Huffman, 1998).Figure 19 presents an example of the incident and scattered fields by a sphere.The second part of the method is an iterative algorithm to calculate scattered field for N spheres from one sphere (figure 20).For the first order, we calculate the scattering of the incident field for each sphere.For the second order, the scattered field of first order for one sphere becomes the incident field for the N-1 other spheres.With this new incident field, the scattered field is calculated as at first order and so on, for higher orders.This iterative process stops when it is converged.The total scattered field is the contribution of all spheres and all orders.The material and the size of these N spheres can be different.Moreover, spheres may be put in a random way, without symmetry conditions on spheres positions.For periodic structures, calculation is simplified.Figure 21 plot the transmission coefficient for the previous infinite dielectric PC.The transmission calculation of 8 layers is calculated in 11 minutes with 63 Mo of working memory.By using the principle of KKR-method of solid state physics, MST can also calculate the band structure (Wang and Al, 1993).MST is the fastest method for the finite structures and random structures.
Photonic Crystals -Introduction, Applications and Theory

Fig. 4 .
Fig. 4. Reduction of the calculation domain in PC

Fig. 5 .
Fig. 5. Calculation of E-field in the PC with PBC, ABC and incident wave

Fig. 7 .
Fig. 7. Band structure and transmission coefficient of dielectric PC

Figure 9
Figure9plot the transmission coefficient for the previous dielectric PC.The transmission calculation of 8 layers is calculated in 32 minutes with 58 Mo of working memory.FDTD is faster than FEM but less accurate.

Fig. 10 .
Fig. 10.Electric field amplitude of dielectric PC at the R-point of the Brillouin zone mesh size along x, y and z-axis.The transfer matrix of the unit cell is written from the unit mesh:

ˆ0T
is the transfer matrix of a vacuum layer.The eigenvectors of ˆ0 T define a plane-

Fig. 11 .
Fig. 11.Band structure and transmission coefficient of dielectric PC

Fig. 13 .
Fig. 13.Band structure of the previous dielectric PC -dimensional case, we choose the constant permittivity along z-axis.The vectors K and G are in the xy plane.The vectors of the central equation are indicated on figure 14.

Fig. 14 .
Fig. 14.Vectors definition of the central equation In the central equation, the matrix of the scalar products is simplified: Photonic Crystals -Introduction, Applications and Theory

Fig. 15 .
Fig. 15.Field expansion on backward wave and forward wave

Figa
Fig. 16.Refraction index of one dimensional PC

.
Similarly to the Kronig-Penney model, the above relation is written from the transmission coefficient: Photonic Crystals -Introduction, Applications and Theory equation 57, we get the multiple-scattering equation:

Fig. 18 .
Fig. 18.Transmission coefficient of triangular 2D-PC with a defect the Hankel and Bessel spherical functions and   , lm Y  are scalar spherical wave functions.The incident field and the scattered field are expanded on spherical wave functions

Fig. 19 .
Fig. 19.Pictures of incident and scattered fields by a metallic sphere Photonic Crystals -Introduction, Applications and Theory

Fig. 20 .
Fig. 20.Block diagram of multiple-scattering method, the total scattered field is the sum of the scattered fields for all orders.

Fig. 21 .
Fig. 21.Transmission coefficient of the infinite periodic dielectric PC Outside the cylinders, the wave equation can be split into two equations, one for the incident field E 0 and the other one for the scattered field E s :