Modeling of the angular tolerancing of an effective medium diffractive lens using combined finite difference time domain and radiation spectrum method algorithms

: A new rigorous vector-based design and analysis approach of diffractive lenses is presented. It combines the use of two methods: the Finite-Difference Time-Domain for the study in the near field, and the Radiation Spectrum Method for the propagation in the far field. This approach is proposed to design and optimize effective medium cylindrical diffractive lenses for high efficiency structured light illumination systems. These lenses are realised with binary subwavelength features that cannot be designed using the standard scalar theory. Furthermore, because of their finite and high frequencies characteristics, such devices prevent the use of coupled wave theory. The proposed approach is presented to determine the angular tolerance in the cases of binary subwavelength cylindrical lenses by calculating the diffraction efficiency as a function of the incidence angle.


Introduction
Subwavelength diffractive optical elements (SWDOE) refer to diffractive optical elements (DOE) that have features sizes smaller than the wavelength of the illumination light. These elements have been studied for several years [1]. They have a growing interest because they can achieve high diffraction efficiency by using only two levels for fabrication [1]. Their "superwavelength" counterparts are diffractive elements that have lateral features larger than the wavelength of illumination. To validate our approach that we propose, we present here the simulation of a subwavelength diffractive cylindrical lens. We demonstrate its ability to focus an incident beam light as a function of the incidence angle.
The diffraction efficiency of a digital scalar DOE is related to the number of phase levels. For example, a 2-levels element, designed by using the scalar theory, can achieve a diffraction efficiency of 40.5% in both first orders. 4-levels diffractive components can reach 81% [2]. The phase quantification which fits the ideal phase function is theoretically better if more levels are fabricated. Therefore, a higher diffraction efficiency becomes possible for a multilevel DOE. This kind of multilevel diffractive element is usually fabricated by multiexposures photolithography technique and requires successive expositions using a high optical precision alignment system. For instance, a 2 N -levels element involves N masking steps. Due to both these successive steps and the reducing of the size of the critical feature sizes of the successive masks, their lateral alignments become more and more critical, since a small misalignment between masks seriously affect the quality of the multilevel DOE [3].
A great interest of binary SWDOE is accounted for their fabrication method. Indeed, one single step of photolithography is required, avoiding the alignment process between masks while maintaining very good diffractive efficiencies [4][5][6][7][8][9]. Prism-like gratings [4], antireflective surfaces [5], computer generated holograms [6], array illuminator [7] or diffractive lenses [8][9][10] designed by effective medium theory have been introduced in the past years. These elements are made up of subwavelength structures that define an effective medium. By varying the lateral feature size, the index of refraction of the effective medium, defined by the subwavelength structures, varies and creates a SWDOE corresponding to a multilevel element. In this case, light "sees" a medium containing an effective modulation of the refractive index. This modulation corresponds to a zeroth order grating which has a high efficiency in the zeroth order, and changes the phase of the incident light [11]. The fill factors of the structures, composing the effective medium, induce a change of the phase of the incoming light. One dimensional [1,9] and two dimensional [4,6,12] surface modulations type have been successfully proposed in the literature in order to create different fill factors and, therefore, different phase modulations.
This paper focuses on a method that proceeds with one dimensional subwavelength DOEs simulations. One of their interests is to avoid the optical alignments of the successive photolithographic masks, while keeping high diffraction efficiency. Our novel approach to simulate subwavelength elements is based on a rigorous propagation of the light through a dielectric structure, using what we called a combined "Finite Difference Time Domain" -"Radiation Spectrum Method" (FDTD -RSM method) algorithms [13,14]. This method enables us to describe the vectorial behavior of light over a long distance. In particular, a f/6 subwavelength diffractive lens has been designed by pulse width modulation and has been studied by the FDTD-RSM method that we propose.

Proposed combined FDTD-RSM numerical method principles
In our proposed approach, we use the FDTD method to calculate the propagation of light through the dielectric structure in the near field and the RSM to simulate the propagation of the light through free space preserving the vectorial aspects of the electromagnetic fields.
The Finite Difference Time Domain is a very well known vector-based method used to simulate the propagation of light by discretising Maxwell's equation on a grid in space and time [13]. To achieve our calculations, we use a FDTD solver called Meep, which is an open software package with subpixel smoothing for increased accuracy, developed at Massachusetts Institute of Technology [15]. Since the studied geometry has only one dimension, we choose to use one longitudinal coordinate for field propagation and one transverse coordinate for the definition of the profile. An advantage of the FDTD simulation is that a large range of different geometries can be simulated, if the lattice is correctly chosen. However, this method, well suited for near field calculations, is very time-consuming. The need to calculate every component of the complex electric and magnetic fields for each space and time grid point limits its use to local scale.
Since the aim is to simulate the propagation of light from the front of the lens to its focal point, which is localized far from the center of the lens, a modal method called Radiation Spectrum Method (RSM) [14] is used. The RSM is a modal Beam Propagation Method (BPM). It is based on the rigorous analysis of the radiation mode spectrum, and has been mainly used in guided optics (e.g. for the reflection calculation at the end of a strongly guiding optical waveguide [16]). The code source enables the determination of both electric and magnetic fields components in the free propagation space resulting from any electromagnetic excitation. Its main advantages concern both the calculation time, the radiation modes sampling conditions and the window edge reflecting problem.
Therefore, in order to make a complete simulation of a SWDOE, jointly concerning near and far fields, the modeling and simulation problem is divided into two parts with appropriate solutions. The first one uses FDTD to evaluate the electromagnetic fields near the diffractive component. The second one describes the propagation in the free space by the RSM up to the focal point of the lens, what considerably reduces the modeling and simulation computing time.

Binary subwavelength design
One characteristic of a binary subwavelength grating under normal incidence is that the only diffracted electromagnetic wave is the zeroth order [11]. In this case, the wave only observes a phase modification related to the geometrical parameters of the grating. In general cases, analytic methods are usually preferred to evaluate the phase difference introduced by this kind of grating [4,10,11]. The drawback is that this kind of methods is limited to basic layout.
In this paper, we propose a numerical method to evaluate the phase difference introduced by a binary subwavelength grating based on FDTD. Thanks to the possibility to solve periodic boundaries conditions problems by Meep, the propagation of a monochromatic wave is calculated over a single period of the grating. From a practical point of view, we consider a TE polarized plane wave which is generated about 5 wavelengths units (roughly 3 µm) below the grating. The TE polarization stands for an electric vector that is perpendicular to the incident plane. Then, the wave propagation is computed through the grating until the steady state is reached. Finally, the phase difference is estimated for all the points located within distance ranging from 1 to 5 wavelength units away from the grating (between 0.6 and 3 µm), and is calculated by estimating the arguments of the electric complex field.
We were looking to simulate a 2 N phase levels surface spatial modulation in order to design the layout of a subwavelength lens which acts as a 2 N levels digital lens. For a dielectric grating in a transmission configuration, the curve giving the phase difference as a function of the fill factor presents a strictly positive slope [11]. Then, the maximum phase difference is reached for an aspect ratio value equal to 1 (i.e. a plane surface corresponding to the height level of the grating). The height of the grating h and the maximum phase difference φ max are linked therefore by Eq. (1): where n is the refractive index of the medium at the wavelength λ. The diffraction efficiency of the zeroth order is defined by the ratio between the electromagnetic energy (calculated by the flux of the Poynting vector) before the grating, and the energy in the zeroth order. For simulations presented in Fig. 1 below, the considered refractive index was 1.46, corresponding to the index of the fused silica at the 632.8 nm designed wavelength that we selected (He-Ne laser). Using the FDTD method, we calculated the main parameters of the zeroth order gratings. An eight phase levels diffractive structure simulated by the effective medium method was targeted, for which the higher phase difference was 7/8 x 2π rad. According to Fig. 1, diffraction efficiency of the zeroth order remains high for gratings with period smaller than 600 nm. In order to obtain the optimum feature size for the subwavelength structure, the phase difference is plotted to link the phase level of the scalar DOE to the fill factor of the subwavelength grating, for a given period. The following graph (Fig. 2) shows the phase difference ∆φ versus the fill factor obtained by FDTD simulations, for several periods. These results are compared to the value computed from relation (2) given by zeroth order Effective Medium Theory (EMT) [11,17] where f is related to the fill factor of the grating and n to the refractive index: For periods smaller than the third of the wavelength, a good agreement between the phase difference computed by the FDTD method and the zeroth order EMT approximation is observed. But the zeroth order EMT is only valid for grating periods smaller than the wavelength.
These results are taken into account to determine eight fill factors corresponding to the required phase differences of an 8-levels diffractive structure. However, our aim is to propose a design that can fit fabrication technologies. The period of 600 nm is only considered, even if the diffraction efficiency of the zeroth order is higher in the case of a 100 nm period grating, what is hard to achieve experimentally. This period leads us to consider the following fill factors and grating sizes, reported in the following Table 1:

Lens design method
The proposed method to design and simulate the behavior of SWDOE is explained in detail in this section. Our aim is to design a diffractive lens operating with subwavelength structures in order to use the coding scheme of the phase.

Diffractive lens layout design
The classical formulae for plano-convex lens, linking the radius of the curvature of the dioptre R, the refractive index n and the focal length f l , is used (Eq. (3)): ( ) The profile of the superwavelength diffractive lens is estimated by "cutting" the phase function into 2π slices with a straight forward design. The resulting profile defines a simple Fresnel lens. Once the 2π modulation has been reached, the resulting profile may be quantified into N phase levels to obtain a N-levels diffractive lens. In fact, the diffractive lens introduces a phase modulation described by a relief surface to the incident wave. As we observed in the previous section, zeroth order gratings can be considered as DOEs which modulate the incident phase. Our aim is now to change the relief modulation into a new structure which contains only binary subwavelength gratings. Figure 3 summarizes the process to design the subwavelength lens.

Electromagnetic wave propagation simulation validating the proposed approach
In order to simulate the electromagnetic behavior of the lens, the combined FDTD -RSM approach that we propose is employed, taking into account all the components of the electric and magnetic fields. Indeed, scalar simulation is not accurate enough because of the subwavelength characteristics of the element.
First, the near field simulation is done with Meep, in order to calculate the propagation of light over a few micrometers after the dielectric structure of the lens. To ensure a sufficient resolution, several samplings are tested and the results compared to each other. Once this propagation is calculated, the components of the fields after the dielectric structures are retrieved in order to become the source of the RSM program.
To calculate the propagation of light through the free space, we use the RSM method. From any electromagnetic source, the guided and radiation modes of a dielectric medium are estimated using this algorithm. In our case, only the radiation modes are considered because the propagation is done through the free space.
The diffraction efficiency η of the lens is then estimated as the flux of the Poynting vector through a surface that corresponds to the diameter of the Airy disk at the focal point. For a lens diameter D, a focal length f, a focal point located at the y f coordinate, at a working wavelength λ and at final step of the simulation t f , the diffraction efficiency is calculated by Eq. (4): , 0 , , .
where × denotes the cross product of the electric E and magnetic H vector fields, x u the u th spatial discrete coordinate x, at the focal plane, and x v the v th spatial discrete coordinate x, at the incident plane, ∆x FDTD and ∆x RSM stand respectively for the sampling steps for FDTD and RSM computation. Figure 4 shows the principle of the "combined" simulation. Left on the Fig. 4, we present amplitude of electric field E inside the dielectric structure and the recovery plane used for the RSM simulation. The Poynting vector is then carried out, step by step, with the RSM to analyse the energy propagation. A detail in the focal plane is observed on the Fig.  4. The diameter of the lens is here 1 mm with a 6 mm focal length. Computation time length of this simulation is approximately 2 hours for FDTD and several minutes for the RSM on Dual Core Pentium(R) D CPU 3.40 GHz Computer model. Fig. 4. Near electric field propagation modeling in an inhomogeneous medium by the FDTD method on the left (the subwavelength lens is located inside the dotted black line) and far field Poynting vector simulation by RSM on the right. The plane A corresponds to the incident plane for the FDTD and the plane B to the recovery plane for the FDTD and also the incident plane for RSM.

Angular tolerance analysis in a specific case
In this section, we present the simulation of a 6 mm focal length focussing subwavelength lenses, described in the previous section. In order to simulate rigorously the behaviour of such lenses, spatial, temporal and spectral resolutions of each simulation steps have to be small enough. Nevertheless, in order to optimize the calculation time, these values should not be too low. Different samplings for the FDTD method were tested to ensure the convergence of simulations for a computational cell size of 1 mm by 5 µm with perfectly matched layer boundaries condition. The obtained fields were compared in order to find the minimum correct sampling. A 20 nm resolution step was found as an enough value to solve properly the propagation of electromagnetic waves through the dielectric subwavelength lenses for the FDTD step, because the electromagnetic fields values were the same for smaller resolutions (mean square error less than 10% for each point). The RSM method has been implemented using a sampling step of 1.3 µm toward Y direction, what defined 800 points in the spatial domain, and therefore, 800 modes in the spectral domain so as to respect the Shannon's sampling theorem. The sampling in the X direction has been 100 µm and 150 reconstructions have been calculated to draw the Poynting vector field around the focal point. In this case, this value is high enough because the total energy is mostly preserved during the propagation, excluding losses due to the divergence out of the calculation window (more than 90% of the incident energy is conserved). Two different implementations of subwavelength lenses were designed and their optical tolerancings were studied. The first one uses 8 phase levels-simulated structures (considering the phase function previously presented in Fig. 2). The second one is a 4 levels simulated lens for which another phase function was calculated, because the maximum desired phase difference was not the same (3π/2 vs 7π/4). As a consequence, the high level value of the gratings composing this 4 levels simulated lens dropped from 1203.5 nm to 1031.6 nm.
In order to determine the angular tolerance of the proposed lenses, the diffraction efficiency was considered as a function of the incident angle of illumination. Figure 5 presents the comparison results:  Figure 5 confirms that the diffraction efficiency of the 8 levels simulated lens is higher than the 4 levels simulated one. Thanks to the scalar theory, it is known that diffraction efficiency is a function of the number of levels. This relation is also verified, because the subwavelength structures are defined thanks to the complete classical diffractive profiles.
The diffraction efficiencies follow the same curve for the angular tolerance whereas the number of phase levels has changed. 65% and higher diffraction efficiencies are obtained for angles smaller than 25°. For angles higher than 40°, diffraction efficiency drop to less than 30%. As the angle increases, no significant differences are observed between the 8 levels and 4 levels cases. So increasing the number of levels is not a rewarding strategy, because the phase sampling is not improved in the case of subwavelength structures. Indeed, at worst, a 8 levels simulated lens requires 1.2 µm high gratings whose line width can be 54 nm thick (see Table 1) although a 4 level simulated lens needs only 1 µm with a 150 nm line thick. That means an aspect ratio of about 22 in the first case versus about 7 in the second case. Therefore, 4 levels seem to be a good compromise between technological constraints and optical performances with angular tolerancing.

Conclusion
SWDOEs have a great interest since they offer the possibility to reduce the fabrication process steps to one single photolithographic step, avoiding optical alignment issues. Thanks to the strong improvement of nanotechnology processes, SWDOEs could be used into micro-electro-opto-mechanical-systems (MEOMS). Therefore, specific methodologies are necessary to design on demand optic functions based on their specificities. This paper presents a new rigorous vector analysis approach fitting these requirements. The proposed algorithm enables to simulate the propagation of light over subwavelength, aperiodic and finite optical elements. The Finite Difference Time Domain algorithm is used to simulate the propagation in the near field condition, and the Radiation Spectrum Method enables the computation of the propagation in the far field. This approach has been used to study a subwavelength diffractive cylindrical lens design, implemented to simulate 4 and 8 levels of phase.
Nevertheless, because 8-levels simulated lens is more difficult to fabricate (because of a higher aspect ratio of the subwavelength gratings for which the width of the line can be smaller), 4-levels simulated lens seems to present the best compromise between good diffraction efficiency and technological constraints.