Generalization and modularization of two-dimensional adaptive coordinate transformations for the Fourier modal method

The Fourier modal method (FMM) has advanced greatly by using adaptive coordinates and adaptive spatial resolution. The convergence characteristics were shown to be improved significantly, a construction principle for suitable meshes was demonstrated and a guideline for the optimal choice of the coordinate transformation parameters was found. However, the construction guidelines published so far rely on a certain restriction that is overcome with the formulation presented in this paper. Moreover, a modularization principle is formulated that significantly eases the construction of coordinate transformations in unit cells with reappearing shapes and complex sub-structures. © 2014 Optical Society of America OCIS codes: (050.1970) Diffractive optics; (160.5298) Photonic crystals; (160.3918) Metamaterials; (050.1755) Computational electromagnetic methods. References and links 1. K. Busch, G. von Freymann, S. Linden, S. F. Mingaleev, L. Tkeshelashvili, and M. Wegener, “Periodic nanostructures for photonics,” Phys. Rep. 444, 101–202 (2007). 2. L. Li, “Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings,” J. Opt. Soc. Am. A 13, 1024–1034 (1996). 3. P. Lalanne and G. M. Morris, “Highly improved convergence of the coupled-wave method for TM polarization,” J. Opt. Soc. Am. A 13, 779–784 (1996). 4. G. Granet and B. Guizal, “Efficient implementation of the coupled-wave method for metallic lamellar gratings in TM polarization,” J. Opt. Soc. Am. A 13, 1019–1023 (1996). 5. L. Li, “Use of Fourier series in the analysis of discontinuous periodic structures,” J. Opt. Soc. Am. A 13, 1870– 1876 (1996). 6. G. Granet, “Reformulation of the lamellar grating problem through the concept of adaptive spatial resolution,” J. Opt. Soc. Am. A 16, 2510–2516 (1999). 7. G. Granet and J.-P. Plumey, “Parametric formulation of the Fourier modal method for crossed surface-relief gratings,” J. Opt. A 4, S145–S149 (2002). 8. T. Vallius and M. Honkanen, “Reformulation of the Fourier modal method with adaptive spatial resolution: application to multilevel profiles,” Opt. Express 10, 24–34 (2002). 9. T. Weiss, G. Granet, N. A. Gippius, S. G. Tikhodeev, and H. Giessen, “Matched coordinates and adaptive spatial resolution in the Fourier modal method,” Opt. Express 17, 8051–8061 (2009). 10. S. Essig and K. Busch, “Generation of adaptive coordinates and their use in the Fourier Modal Method,” Opt. Express 18, 23258–23274 (2010). 11. J. Küchenmeister, T. Zebrowski, and K. Busch, “A construction guide to analytically generated meshes for the Fourier Modal Method,” Opt. Express 20, 17319–17347 (2012). 12. J. Küchenmeister, “Three-dimensional adaptive coordinate transformations for the Fourier modal method,” Opt. Express 22, 1342–1349 (2014). #207705 $15.00 USD Received 5 Mar 2014; revised 1 Apr 2014; accepted 2 Apr 2014; published 10 Apr 2014 (C) 2014 OSA 21 April 2014 | Vol. 22, No. 8 | DOI:10.1364/OE.22.009404 | OPTICS EXPRESS 9404 13. L. Li, “Fourier modal method for crossed anisotropic gratings with arbitrary permittivity and permeability tensors,” J. Opt. A 5, 345–355 (2003).


Introduction
Periodic arrangements of nanostructures haven proven to be a rich field of research [1].Alongside the rapid development in experimental techniques, the methods to solve Maxwell's equations numerically have significantly evolved as well.
One of the solvers that has proven to be of great applicability is the Fourier modal method (FMM).It predicts the transmission properties of periodic photonic systems, both dielectric and metallic.The considered systems are periodic in one or two directions (xy-plane) and finite in z-direction.The system is sliced into layers with constant permittivity in z-direction, each of which lead to an eigenvalue problem representing Maxwell's curl equations.This allows expanding the fields into eigenmodes.The layers are later joined together using a scattering matrix algorithm.This algorithm also ensures that the continuity conditions are fulfilled [2].
A major step forward was achieved when the correct Fourier factorization rules [3][4][5] were found.However, metallic systems were still hard to solve.The basic problems of in-plane staircasing for not-grid-aligned structures and most of all the Gibbs phenomenon still had to be tackled.
An elegant way to address these problems was found in coordinate transformations.Two categories were designed: Firstly, adaptive coordinates (AC) transform the permittivity distributions such that the surface of the considered structures becomes grid-aligned.Secondly, the coordinate line density at the structure's surface is increased with adaptive spatial resolution (ASR), another kind of coordinate transformation.When both are combined, a great variety of structures can be reliably treated with the FMM [6][7][8].In the last years, different concepts were found for the construction of the corresponding meshes in two [9][10][11] and even in three dimensions [12].
In particular, [11] introduced a recipe how to construct meshes for a great variety of structures in two dimensions.However, a certain restriction in the formalism does not allow arbitrary shapes of the mesh.The main scope of this paper is to eliminate this restriction and provide a short mathematical form to construct arbitrary meshes.In order to do so, we first discuss the basics of Maxwell's equations in generalized coordinates in section 2.Then, we proceed to the generalization of the construction procedure in section 3. Here, the restriction in [11] is described and solved for a simple test case.Thereafter, the formulation for arbitrary shapes is presented.Finally, section 4 presents a modularity concept for mesh construction.It constitutes an option for unit cells with repeating shapes within the unit cell or shapes with complex substructures.For those systems, a way to incorporate previously known meshes is outlined that shortens the overall mesh procedure.
The reader that is only interested in the formulation of the fully generalized formulation may skip ahead to Eqs. (15) to (20).

Covariant formulation of the Fourier modal method with generalized coordinates
This section discusses the impact of coordinate transformations on Maxwell's equations and the consequences for the Fourier modal method.Since this is not the first publication on this topic, I keep the section short and use the notation already used in [9][10][11][12][13].The type of system considered in this paper in displayed in Fig. 1, one layer of periodic nanostructures illuminated from above.
In the following, we use a curvilinear coordinate system Ox 1 x 2 x 3 and a Cartesian coordinate system Ox 1 x 2 x 3 .The coordinate transformations we consider in this paper are two-dimensional Fig. 1.Systems like a layer of rotated gold crosses cannot be simulated reliably with classical FMM.For that, adaptive coordinates and adaptive spatial resolution have shown to be a powerful tool to improve the representation of the structure.The efficient mesh construction for this type of system is the scope of this paper.and have the form (1) Maxwell's equations in covariant form read [13] Here, E σ and H σ are covariant components of the electric and magnetic field and ξ denotes the Levi-Civita symbol.Einstein's sum convention is used (repeated indices are implicitly summed over) and the Greek indices run from 1 to 3. The vacuum wave number is denoted k 0 = ω/c with the frequency ω and the speed of light c.The metric tensor g ρσ reads where g (as used in Eqs. ( 4) and ( 5)) is the reciprocal of its determinant.Applying a coordinate transformation leads to the so-called effective permittivity Here, ε τχ is the permittivity tensor in the Cartesian system.The effective permittivity is the matrix that is Fourier transformed in the FMM.Therefore, the aim of coordinate transformations is to obtain an effective permittivity that is ideally suited for the Fourier transformation.
Particularly, this means that the surfaces are grid-aligned in the effective permittivity and the coordinate line density at the surface of the structure is increased (see [11] for details and illustration).

Generalization of two-dimensional mesh construction
In this section, I first briefly sketch the construction principle from [11] and describe the restriction implied in this scheme.Secondly, I demonstrate how this can be overcome using a simple example.Finally, I generalize the construction principle to arbitrary functions.

Description of the problem
The mesh construction principle from [11] is sketched in Fig. 2. The example shown is a unit cell with a crescent-shaped structure.The rest of the mesh in panel (c) is a linear transition between the characteristic coordinate lines and the unit cell edge.This mesh construction scheme was presented in [11] where the mathematical formulation contains the constraint that the characteristic points P, Q, R, S need to be mapped onto themselves.This constraint will be overcome in the following.The figures are taken from [11] and are slightly modified.
To obtain a mesh, the so-called "characteristic points" were chosen on the surface of the structure, see P,Q,R,S in Fig. 2(a).The characteristic points in turn defined characteristic coordinate lines which were mapped in a way that they match the surface of the structure, cf.blue and red lines in Figs.2(a) and 2(b).To obtain the total mesh, the mapping in between the characteristic coordinate lines was performed as a linear transition between the edges of the unit cell and the characteristic coordinate lines.The mesh resulting from this procedure is depicted in Fig. 2(c).
In the mathematical formalism in [11] an implicit restriction was posed: the characteristic points had to be mapped onto themselves.As it turns out, there are cases in which this restriction limits an efficient mesh construction.An example is the structure shown in Fig. 1.Here, the simplest way to construct a mesh is depicted in Fig. 3.
Unfortunately, this way of constructing the mesh relies on characteristic points that are not mapped onto themselves.An important goal of this paper is to formulate a generalization of the formalism in [11] that allows an arbitrary choice of the characteristic points' mapping and an arbitrary choice of the characteristic coordinate lines' mapping.The following section deals with the formulation for the mesh derived from the sketches in Fig. 3, which is a special case due to the fact that the functions connecting the characteristic points are linear.

Solution for piecewise linear mappings
For the discussion of the mapping, we recall the LT (c, c, d, d, x)-function introduced in Eq. ( 13) in [11].It defines a straight line through the points (c, c) and  show how to obtain a mesh for the structure in Fig. 1 where the characteristic points are not mapped onto themselves, i.e., P 0 = P ϕ .Using the formulation from Section 3.2, the mesh in panel (c) can be constructed.The sketches are taken from [12] where the construction principle presented in Section 3.2 was already applied.
We only discuss the mapping for zone 1 in Figs. 3 and 4(a), the rest follow in analogy.The mapping in the first zone reads with with Let us interpret Eqs. ( 9) to (12): in Eq. ( 9) the unit cell edge is mapped onto itself (0,0) and the P 0,x 1 coordinate line, blue in Fig. 3(a), is mapped using the help function H 1 (x 2 ).Something similar happens in Eq. (11) where also the unit cell edge is mapped onto itself (0,0) and the P 0,x 2 coordinate line, red in Fig. 3(a), is mapped using the help function H 2 (x 1 ).We observe that the mesh construction is conceptually new in the help functions H i which are visualized in Fig. 4(b).The help functions are not the functions that directly connect (P 0,x 1 , 0) and (0, P 0,x 2 ) with P ϕ which is the essential change in the mapping procedure compared to Section 3.1 and [11].The point is that we cannot map each coordinate line independently from the other coordinate line any more since we have to make sure that P 0 is mapped onto P ϕ .This means, that the x 1and x 2 -mappings need to be designed such that they yield the desired results when they are combined.This is accomplished by the form of the above mapping.We also note that the LT functions in Eqs. ( 9) and ( 11) describe the linear transition between the characteristic coordinates lines.This means that these LT functions define the coordinate line density in each zone.The help functions H i , on the other hand, are linear since the characteristic coordinate lines are mapped on piecewise linear functions, see Eqs. ( 10) and ( 12) and Fig. 3(b).The next goal is to generalize this approach to allow any shape at the edge of the mapped zone.Fig. 4. Construction of the help functions.To map the zone 1 from Fig. 3 (shown in panel (a)) correctly, the help functions H 1 and H 2 need to be chosen such that they both combined yield the correct result.Especially, we have to make sure that the mapping procedure maps P 0 on P ϕ .

General construction scheme for arbitrary help functions
In order to cover arbitrary functions like g and h in Fig. 5(a), we have to change the help functions H 1 and H 2 .The mappings in Eqs. ( 9) and ( 11) stay like they are (which is in line with the scheme discussed in [11]) and we only vary the help functions to successfully construct the overall mesh.g h Fig. 5. Similar to Fig. 4 the characteristic point is not mapped onto itself.In addition, the mapped coordinate lines are not linear any more.To construct the suitable mapping, the help functions need to be varied as discussed below.The resulting mesh for this zone is displayed in panel (b).
As stated above, we need to make sure that P 0 is mapped on P ϕ and that the mappings in both directions are combined correctly.The formulation for this mapping can be condensed to the following form: which is plugged in Eqs. ( 9) and ( 11) to obtain the full mapping in zone 1 .The intermediate step is supposed to emphasize that we cannot just use the functions h and g.Instead, the stretched or compressed coordinates x 2 h and x 1 g make sure to map P 0 on P ϕ .This formulation is rather elegant since it does not matter where P ϕ is in relation to P 0 -the formulation holds regardless of whether the zone is e.g.stretched in one direction and compressed in the other.The resulting mesh constructed with a quadratic g function and a cubic h function is depicted in Fig. 5(b).
After discussing the previous special cases, I now present the most general formulation for a zone where all four vertices and the connecting functions are arbitrarily chosen.This is depicted in Fig. 6.In Fig. 6(a) we see the rectangle that is supposed to be mapped.The new vertices and the connecting functions are depicted in Fig. 6(b).Once the formulation for this situation is known, any mesh can be easily constructed.We start by giving the transformation similarly to Eqs. ( 9) and (11).It reads Equation ( 15) refers to the vertical coordinate lines and expresses that the mapping is a linear transition between the two blue coordinate lines in Fig. 6(b).For that, the coordinate line P 1,x 1 is mapped on e and the S 1,x 1 coordinate line is mapped on g.Similarly, Equation (16) states that the P 1,x 2 is mapped according to the function f , the coordinate line Q 1,x 2 is mapped according to g, cf.vertical, red lines in Figs.6(a) and 6(b), and again a linear transition is performed between them.As stated above, we have to take into consideration that the vertices are not mapped onto themselves, which results in a transformed coordinate.This is taken care of by the following formulation Figure 6(c) shows an example of a mesh constructed with this formulation.Here, the function f connecting the points P 2 and S 2 was chosen to be a polynomial function of fourth order.Once this function was known, Eq. ( 18) could be written explicitly and plugged into Eq.(15).
All in all, Eqs. ( 15) to (20) form an easy recipe for any user of adaptive coordinate transformations.Only the vertices of the desired mappings and the connecting functions need to be known and plugged into Eqs.( 15) to (20) to obtain the desired mapping.

Modularized construction of two-dimensional meshes
With the formalism discussed in the previous section, we have generalized the mesh construction to arbitrary shapes.This section, on the other hand, is devoted to a practical problem: in many systems of interest, certain shapes occur several times within the unit cell or the unit cell consists of an arrangement of previously meshed structures.For these situations, it makes sense to set up certain modules that can be reused several times in one mesh.The basic modularity concept is illustrated in Figs.7 and 8. Fig. 7.The modularization principle can be used when a mesh for a structure with complex sub-structures needs to be found.To obtain it, the characteristic points and coordinate lines in panel (a) are chosen such that the area for the sub-structure is mapped onto itself, see area A in panel (b).The resulting preliminary mesh is depicted in panel (c).In the second step, the sub-structures are entered in this region, see Fig. 8.
Here, we define characteristic lines that are quite similar to the ones in Fig. 7 in [11], i.e., a mesh for a large circular structure with a sub-structure in the middle.The difference to the mesh in Fig. 7 in [11] is an altered behavior in the center: instead of meshing a specific substructure like another circle, we intentionally map the area A onto itself.Thereby, we have a Cartesian area in the middle of our mesh.This is where the idea of modularity comes into play-in a second step we can replace the Cartesian mapping with any of the mappings that we have acquired already (properly scaled, of course).An example is shown in Fig. 8(a), where the Cartesian mapping in the center of the large circle is replaced with a small circle.Alternatively, we could have put in the crescent or any other meshed structure.Moreover, we can assemble even more complex structures-this is depicted in Fig. 8(b), where we entered several more circles.Using modularity and entering properly scaled (and previously known) mappings eases the construction significantly.In contrast, imagine the large number of characteristic points and coordinate lines that would have to be considered for the structure in Fig. 8(b).

Conclusion
In this paper, the two-dimensional adaptive mesh construction for the Fourier modal method was generalized far beyond the scope of [11].The formulation presented here allows for an arbitrary choice of the characteristic points' mapping and of the functions connecting these points.A closed formulation of the principle was presented that allows to quickly reproduce and adapt the formalism.In addition to the generalization of two-dimensional meshing, a principle for modularization was presented that allows to construct meshes that combine previously known shapes and form complex sub-structures without much effort.

Fig. 2 .
Fig.2.The characteristic coordinate lines in panel (a) (red and blue) are mapped onto the surface of the structure, see panel (b).The rest of the mesh in panel (c) is a linear transition between the characteristic coordinate lines and the unit cell edge.This mesh construction scheme was presented in[11] where the mathematical formulation contains the constraint that the characteristic points P, Q, R, S need to be mapped onto themselves.This constraint will be overcome in the following.The figures are taken from[11] and are slightly modified.

Fig. 3 .
Fig.3.This figure demonstrates that a different choice of the characteristic points can be reasonable.Panels (a) and (b) show how to obtain a mesh for the structure in Fig.1where the characteristic points are not mapped onto themselves, i.e., P 0 = P ϕ .Using the formulation from Section 3.2, the mesh in panel (c) can be constructed.The sketches are taken from[12] where the construction principle presented in Section 3.2 was already applied.

Fig. 6 .
Fig. 6.Generalized construction principle.The rectangular area in panel (a) is mapped in an arbitrary fashion in panel (b).The characteristic points are not mapped onto themselves and no further restriction is imposed.The connecting functions e, f , g, h are chosen freely as well.The mesh obtained with the generalized expression presented below is shown in panel (c).

Fig. 8 .
Fig. 8.In the second step for the modularity concept, the area A that was mapped onto itself, see Figs. 7(a) and 7(b), is replaced by a scaled transformation of a previously known structure, see panel (a).This method also allows for a number of sub-structures as shown in panel (b).