Covariant Description of Transformation Optics in Linear and Nonlinear Media

The technique of transformation optics (TO) is an elegant method for the design of electromagnetic media with tailored optical properties. In this paper, we focus on the formal structure of TO theory. By using a complete covariant formalism, we present a general transformation law that holds for arbitrary materials including bianisotropic, magneto-optical, nonlinear and moving media. Due to the principle of general covariance, the formalism is applicable to arbitrary space-time coordinate transformations and automatically accounts for magneto-electric coupling terms. The formalism is demonstrated for the calculation of the second harmonic generation in a twisted TO concentrator.


Introduction
The field of transformation optics (TO) has drawn a lot of scientific interest in the last few years [1][2][3][4][5][6]. By this design methodology, the form-invariance of Maxwell's equations under coordinate transformation is used to tailor the optical properties of an electrodynamic medium. The majorities of TO applications focus on the design of the linear material parameters. In that case, a coordinate transformation is used to engineer the linear constitutive parameters of a medium such that the wave trajectory follows a desired path [7]. An interesting application of this concept is the realization of electromagnetic invisibility cloaks [8][9][10][11][12]-devices in which light is guided around a certain region of space rendering the interior of the region invisible for an external observer. Cloaking devices belong to the most prominent applications in which the TO concept is successfully applied for designing the linear material properties of a medium and have been extensively reported in the literature [13][14][15][16][17][18][19].
In contrast, only little work in the research of TO media addresses the transformation of nonlinear material properties. Media with nonlinear response, however, provide a number of interesting effects including sum-and differencefrequency generation, parametric amplification and oscillation, stimulated scattering, self-phase modulation and self-focusing [20][21][22][23][24][25][26][27][28] and play a key role in modern optical technology [29][30][31]. Consequently, the extension of the TO concept to nonlinear media is expected to offer a variety of new opportunities for engineering optical media and for the construction of novel electromagnetic devices [32]. Apart from that, the TO concept also provides a promising calculation method in modeling complex nonlinear systems whenever it is possible to find a coordinate transformation such that the geometry of the system takes a much simpler form in the transformed space.
A basic prerequisite for a successful integration of nonlinear effects into the TO concept is a general transformation law for linear and nonlinear material parameters under space-time coordinate transformation. A first step in this direction is suggested in [32] where the transformation of certain classes of nonlinear materials under purely spatial transformations is studied. The aim of the following paper is to generalize this approach and to provide a rigorous theoretical framework for arbitrary nonlinear materials and for both temporal and spatial coordinate transformations. The formalism is presented in a manifestly covariant form that allows the simultaneous treatment of electric, magnetic and magneto-electric cross-coupling effects and encompasses all types of space-time transformations with a particular reference to moving media [33][34][35].
The paper is organized as follows: In the first part, we introduce the used tensor notation and explain how the defined tensors are related to the common linear and nonlinear material parameters such as the permittivity, permeability or the quadratic electro-optic coefficients. Subsequently, we exploit the fundamental principle of relativity to derive a general transformation law for nonlinear media under continuous space-time coordinate transformations. In this context, a special emphasis is given to the transformation of nonlinear constitutive parameters in moving media and the class of time-independent, spatial transformations which play a central role for the design of TO devices. In the final part of the paper, we illustrate and explain the derived expressions by means of an explicit calculation example. By deriving the second harmonic generation in a nonlinear, twisted field concentrator, we show that an appropriate coordinate transformation can provide a significant alleviation in the formal treatment of a complex nonlinear problem and, thus, allows a convenient calculation of an otherwise sophisticated process.

The covariant material equation
In order to establish a common basis for the following discussion and to introduce the used notation, we start with a short review of the covariant description of the electrodynamic theory. In this paper, all quantities are expressed in Gaussian units. In this case, the Maxwell equations in matter take the form: where E, H, D, B, ρ, j and c denote the electric and magnetic field, the electric and magnetic flux density, the charge, the current density and the speed of light in vacuum, respectively. By using the four-current (j ν ) = (cρ, j), the antisymmetric field strength tensor and the antisymmetric displacement tensor the Maxwell equations can be covariantly expressed as: where we use the Einstein summation convention † and the contravariant derivation (∂ µ ) = ( 1 c ∂ t , −∇). In the covariant notation, (5) contains the two homogeneous Maxwell equations whereas (6) contains the two inhomogeneous Maxwell equations. Note that the signs used in the tensor notation above depend on the convention used for the metric tensor. In this paper, we use the metric tensor given by η µν = diag(1, −1, −1, −1). The relation between the tensors F µν and D µν is given by with the polarization-magnetization tensor Tensor equation (7) is equivalent to the common material equations: In many cases, the quantitative relation between the polarization-magnetization tensor P µν and the field strength tensor F µν can be expressed in a power series according † Throughout this paper, Greek indices run from 0 to 3 while Latin indices run from 1 to 3. to: Note that in this notation, the index pair µν is fixed while the other indices are summed over. From the antisymmetry of the tensors F µν and P µν and the commutativity of the products F α i β i F α j β j , it follows that the coefficients χ α 1 β 1 ···αnβn µν are antisymmetric under exchange of µ ↔ ν, symmetric under exchange of two pairs α i β i ↔ α j β j and antisymmetric under exchange of two indices α i ↔ β i within one pair α i β i . For the quadratic term, this means for example: It is obvious that additional, inherent symmetries (such as the Kleinman symmetry or spatial symmetries given by the point symmetry class of the medium) further reduce the number of independent coefficients. However, in order to provide a general description, we only consider the basic symmetries given in (11). Inserting (10) into the material equation (7) yields: where in the last line, we performed a re-definition of the linear coefficient χ α 1 β 1 µν to include the free-space space contribution (see appendix).

The linear term
To become familiar with the tensor notation, it is instructive to first consider only the linear term on the right-hand side of (12): By using the symmetry properties of χ σκ µν , it follows for the components of the electric displacement field D: where ǫ j i denotes the permittivity tensor and ξ j i is the tensor of the electro-magnetic coupling (remember, Latin indices run from 1 to 3). The last line follows from the fact that F 0j addresses the components of the E-field while F kl addresses the components of the B-field. In a similar way, one finds for the components of the H-field: with the totally antisymmetric Levi-civita tensor g mn i (with g 23 1 = 1), the permeability µ j i and the magneto-electric coupling tensor ζ j i . In summary, the linear tensor equation D µν = 4πχ σκ µν F σκ can be expressed in the common three-formalism as: This is the general equation of a linear bianisotropic medium. The relations between the parameters ǫ, µ, ξ and ζ and the four-tensor χ σκ µν can be derived by comparing the coefficients occurring in (14) and (15). As shown in the appendix, one finds:

Nonlinear terms
Next, we consider quadratic contributions to the polarization-magnetization tensor P µν , i.e. contributions that depend on the product of two components of the field strength tensor. These are given by the second summand in (10) and have the form: For simplicity, we restrict to the components of the second-order electric polarization P 0i (a similar derivation applies for the magnetization). As for the linear term, we can use the symmetry properties of χ σκαβ µν to split the summation into three terms: Pockels effect, multi-wave mixing where we have introduced the nonlinear material coefficients a jk i , b jk i , c jk i . In other words, the expression P (2) µν = χ σκαβ µν F σκ F αβ describes all second-order electric, magnetic and magneto-electric cross-coupling effects in a single equation. Similarly, one finds for the nonlinear polarization of the third order: and so on. Note: while the quadratic contributions in (18) contain 6 + 9 + 6 = 21 independent coefficients, the cubic contributions in (19) contain already 10 + 18 + 18 + 10 = 56 independent coefficients (see appendix).

Coordinate transformation
The particular advantage of the covariant formulation is its form-invariance under coordinate transformations where (x α ) = (ct, r) is the coordinate vector of the space time. By A α ′ α = ∂x α ′ /∂x α and |A| := det(A α α ′ ) we denote the Jacobian matrix and its determinant. To clarify the following discussion, we introduce the convention that, if a product of A α ′ α occurs in a transformation formula, the kernel symbol A is written only once, e.g.
The electromagnetic fields F µν and D µν in the primed and unprimed coordinate systems are related by [33,Chap. 3.2]: The Maxwell equations (5) and (6) have the same form in the primed coordinate system as in the unprimed system due to their natural form-invariance [33,Chap. 3.2], that is: In order to achieve form-invariance also for the constitutive relation (12), the material coefficients must transform as: The proof of (24) is not difficult. To see this, we multiply both sides of (12) by |A| −1 A µ ν µ ′ ν ′ , insert in the right-hand side the identity (see [33], Chap. 1.4) and replace the unprimed quantities with the corresponding primed ones. This yields: Thus, the constitutive equation transforms covariantly as required. Relation (24) represents the general transformation law for linear and nonlinear materials parameters (including complicated magneto-electric coupling terms) and is equally valid for arbitrary space-time coordinate transformations.

Moving media
A typical physical situation where coordinate transformations play a role occurs when a medium moves relative to the observer. In such cases, the transformation law allows the calculation of the material parameters in the reference frame of the observer if the material parameters are known in the rest frame of the medium. Consider, for a example, a linear, homogeneous medium that is isotropic in its rest frame with the constitutive equations: According to (17), the covariant material tensor χ σκ µν of an isotropic medium is related to the material parameters ǫ and µ by the equations: (28) (modulo interchange of indices, e.g. ǫ = −8π χ 0i i0 , etc.) while χ σκ µν is zero otherwise. By means of the transformation rule of (24), we can now calculate the material equation in a reference frame in which the medium moves uniformly, say along the x-direction, with constant velocity v. The corresponding coordinate transformation is given by with β = v/c and γ = (1 − β 2 ) −1/2 . The inverse of the Jacobi matrix A α α ′ is obtained by replacing β by −β (corresponding to a velocity reversal). The material equation in the frame in which the medium is moving is obtained by applying the transformation law (24) to the material equation (13) (expressed in primed coordinates), that is After a straightforward summing over the indices and by using the relation between the field tensors D µ ′ ν ′ and F σ ′ κ ′ and the E ′ , B ′ , D ′ and H ′ fields (see Eqs. (3) and (4)) as well as the relation between the material tensor χ σκ µν and ǫ and µ given by (28), this can be rewritten in the more intuitive form: in agreement with the findings in [36]. Since these equations are form-invariant under spatial rotation of the coordinate system, they are equally valid for arbitrary directions of the medium's velocity v. Obviously, a linear medium that is isotropic at rest becomes bianisotropic if it is moved relative to the observer [35][36][37][38].
The mixing of the electric and magnetic response in moving media is a general effect which can also be observed for the nonlinear material coefficients. Consider, for example, the quadratic term of the material tensor given by χ σκαβ µν . In the frame in which the medium is moving the tensor transforms as As before, we can evaluate the sum, collect corresponding terms and identify the new nonlinear coefficients. However, in this generality the final expressions get very complicated and unintuitive. Therefore, in order to illustrate only the principle of the mixing of higher-order susceptibilities in moving media, we restrict to the second-order contribution to the electric polarization P (2) i = χ σκαβ 0i F σκ F αβ . As described in (18), the relevant components of the material tensor can be further split in terms with χ 0n0l 0i , χ 0nkl 0i and χ mnkl 0i (modulo interchange of indices) to separately account for contributions that are proportional to E i E j , E i B j and B i B j , respectively. If we now apply the transformation (32) separately to these three coefficients, the result can be expressed in a compact matrix form: Obviously, the time-dependence A i ′ 0 = c −1 ∂x i ′ /∂t = 0 due to the movement of the medium implies nonzero off-diagonal elements which induce a mixing of nonlinear material coefficients of electric and magnetic nature. This means, for example, that a Pockels medium at rest (χ 0n0l 0i = 0) can display a Faraday effect (χ 0 ′ n ′ k ′ l ′ 0 ′ i ′ = 0) if it is moved relative to the observer and vice versa.

Time-independent transformations
We now focus on the special case of time-independent, spatial transformations, i.e. transformations with which are of particular interest for the design of TO devices with tailored optical properties. In this case, we have Consequently, if the term A α ′ i T α ′ β ′ occurs in an expression, we can replace the fourindex α ′ by the three-index i ′ (which runs from 1 to 3) since A 0 ′ i = 0. For example, the permittivity transforms under these conditions as Accordingly, the permeability and magneto-electric coupling terms in (17) transform as: Note that for purely spatial transformations there is no magneto-electric cross-mixing. By the general transformation law (24), we can now calculate the transformation behavior of the nonlinear material coefficients under spatial coordinate transformations. For instance, we exemplarily calculate the second order susceptibility of the electric polarization a jk i = 4χ 0j0k 0i which describes nonlinear optical effects such as secondharmonic generation or three-wave-mixing. With the help of (24), the general transformation of the second order tensor χ σκαβ µν is and for the special case of a spatial coordinate transformation (i.e. with the conditions given in (35)), the second order susceptibility transforms as where we applied similar calculation steps as in (36). This shows that as soon as the relation between the material coefficient of interest (linear or nonlinear) and the general covariant tensor χ α 1 β 1 ···αnβn µν is found, the transformation law immediately follows from (24).

Twisted, nonlinear field concentrator
In the following, we demonstrate that nonlinear problems in complex media with sophisticated linear and nonlinear optical properties can take a much simpler form in an appropriately transformed space.
To provide an illustrative example of this calculation method, we consider a nonlinear, inhomogeneous material that is illuminated by a strong laser field. We assume that the polarization of the fundamental wave and the nonlinearity of the material match the condition for second harmonic generation (SHG) and, as a goal, we want to calculate the spatial field distribution of the SHG wave inside the medium. If the material is highly inhomogeneous, both the fundamental and the SHG wave follow a complicated, distorted trajectory through the medium which generally hampers a numerical calculation and reliable prediction of the SHG progress. However, if it is possible to find a coordinate transformation such that the wave propagates uniformly along straight lines in the new coordinate system, the wave equation can be readily solved by an ordinary integration along the field lines. A subsequent back-transformation then yields the SHG field in the physical space. In the following, we demonstrate and explain this calculation method for a specific example.
As a hypothetic, inhomogeneous material, we consider the special case of a TO medium, i.e. an artificial material whose permeability and permittivity tensors were obtained by applying a coordinate transformation to a homogeneous, isotropic space. We suppose that the optical properties of the homogeneous space are similar to that of vacuum with a permittivity and permeability equal (or close) to unity. In the following, the coordinates of the uniform, isotropic space are indicated by unprimed indices while the coordinates of the physical space of the inhomogeneous material are indicated by primed indices.
For the transformation between the primed and unprimed system, we consider the following transformation (expressed in cylinder coordinates using r 2 = x 2 + y 2 and φ = arctan(y/x)): As illustrated in Fig. 1, the transformation compresses space of a cylindrical region with radius R 2 into a region with radius R 1 at the cost of an expansion of space between R 1 and R 3 [15]. In addition, the space experiences a radius-dependent twist about the z-axis [39]. Note that the transformation addresses only the x-and y-coordinates while the z-coordinate remains unchanged. For this reason, we can restrict the following discussion to the xy-submanifold of the space-time manifold. The corresponding Jacobi matrix is: And the determinant of A g is: The transformation (40) is expressed in cylindrical coordinates. For the calculation of the SHG, however, it is more advantageous to formulate the transformation in Cartesian coordinates in the form x ′ = x ′ (x, y) and y ′ = y ′ (x, y). This is achieved by applying an intermediate transformation r = x 2 + y 2 and φ = arctan(y/x) in the unprimed system and an inverse intermediate transformation x ′ = r ′ cos φ ′ and y ′ = r ′ sin φ ′ in the primed system. With these expressions, the transformation (40) can be re-expressed in Cartesian coordinates by applying the following three subsequent transformations: According to the chain rule in higher dimension, the Jacobian matrix of a composite function is just the product of the Jacobian matrices of the composed functions, that is As shown in the appendix, the Jacobian determinants of the intermediate transformations are |A f | = 1/r and |A h | = r ′ , respectively. Consequently, the Jacobian determinant of the composite is with the Jacobian determinant |A g | given in (42). Once the Jacobian matrix and its determinant are known, we can easily transform the fields and material parameters from one coordinate system to the other.
We now intend to calculate the SHG wave generated in the concentrator. For this purpose, we suppose that the material used for the construction of the concentrator exhibits a non-vanishing nonlinear susceptibility. To simplify matters, we further assume that this nonlinearity obeys the phase matching condition for frequency doubling if the fundamental wave and the second harmonic wave are both polarized in the z-direction. In our notation, the corresponding tensor component of the nonlinear susceptibility is a ′ : (18)). For the spatial distribution of a ′ , we suppose that the nonlinearity is only present in the inner cylindrical region of the material according to: This could, for example, be realized by doping the center of the concentrator with some nonlinear material.
As fundamental wave we assume an incident monochromatic plane wave that initially propagates along the x-axis while the electric field vector is polarized in the z-direction. In the physical space of the inhomogeneous material (spanned by primed coordinates), the fundamental wave takes the form: where E ′ ω denotes the amplitude of the wave and k ′ ω = (k ′ ω x ′ , k ′ ω y ′ ) is the wave vector in the medium. For the second harmonic wave, we represent the electric field by with the SHG wave vector k ′ 2ω = 2k ′ ω (phase matched case). In the slowly varying amplitude approximation for the second harmonic wave and the undepleted pump approximation (i.e. E ′ ω is constant), the wave equation for the SHG amplitude E ′ 2ω (x ′ , y ′ ) is given by [31]: with κ = −16πiω 2 /c 2 . This is a partial differential equation in two dimensions for the second harmonic field amplitude E ′ 2ω (x ′ , y ′ ) which is very difficult to solve in general. However, the complexity can be significantly reduced if we apply a coordinate transformation to the uniform space (spanned by unprimed coordinates).
According to (22) and the Jacobian matrix of the transformation given by (44) (with suppressed subscript f • g • h), the z-component of the electric field transforms as since t ′ = t and z ′ = z implies A µ ν 0 ′ 3 ′ = δ µ 0 δ ν 3 . Consequently, the electric fields of both the fundamental and SHG wave remain unchanged under coordinate transformation (but the functional dependence on the underlying coordinate system may be different). Since the unprimed system is homogeneous and isotropic with a refractive index equal to one, the electric field of the fundamental wave in the unprimed system has the form of a uniform plane wave propagating in the x-direction: with constant amplitude E ω = E ′ ω and k ω x = ω/c. From this expression, the electric field in the primed coordinate system is readily obtained by expressing x by x = x(x ′ , y ′ ). The real part of the fundamental wave in the primed and unprimed system is plotted in Figs. 2(a) and 2(b), respectively.
Since the fundamental wave propagates straightly along the x-direction in the unprimed system and since we assumed phase matching for the SHG process, it follows that the wave vector of the SHG wave is given by k 2ω x = 2k ω x = 2ω/c and k 2ω y = 0 (on the considered length scale, beam divergence due to diffraction can be neglected). Consequently, the wave equation (49) for the SHG process reduces in the unprimed system to the simple expression which can be immediately integrated to: The remaining unknown is the nonlinearity a = 4χ 0303 03 , i.e. the relevant tensor component for the SHG process in the unprimed system.
According to the transformation law derived in (39), the nonlinearity transforms as: With the determinant |A| given by (45) and (42), the nonlinearity a ′ defined in (46) and the relation between r and r ′ given by (40), we can calculate the nonlinearity in the unprimed system to be: As expected, the nonlinearity in the unprimed system is reduced by a factor of (R 1 /R 2 ) 2 compared to the nonlinearity in the primed system (see (46)) because the cylindrical region in which the nonlinear substance is located experiences a space expansion from radius R 1 to R 2 if we perform a coordinate transformation from the primed to the unprimed system. Now we can (even analytically) evaluate the integral (53) for the amplitude E 2ω (x, y) of the SHG wave. A subsequent multiplication with the propagation phasor yields the SHG wave in the unprimed system in the form: By applying the inverse transformation (i.e. expressing x and y by x = x(x ′ , y ′ ) and y = y(x ′ , y ′ )), we finally obtain the SHG wave in the physical space of the primed coordinates. The real part of the resulting SHG field distribution in the two coordinate systems is plotted in Figs. 2(c) and 2(d), respectively. The proposed twisted nonlinear concentrator is certainly a somewhat constructed example since the exploited uniform space was already presumed in the design of the concentrator. However, the decisive step in the calculation-the straightening of the wave trajectories inside the medium by applying an appropriate coordinate transformation-is in principle always possible in any inhomogeneous media with continuously varying material properties. The proposed technique of simplifying nonlinear processes in inhomogeneous media is therefore not restricted to the presented example, but covers a wide application range.

Conclusion
We proposed a theoretical framework for the incorporation of nonlinear effects within the concept of transformation optics (TO). In this context, we derived a general expression for the calculation of linear and nonlinear electromagnetic material parameters under arbitrary coordinate transformations. The transformation law is formulated in a manifestly covariant form that allows the simultaneous treatment of electric, magnetic and magneto-electric cross-coupling terms and is applicable to both temporal and spatial transformations.
As a first application example, we calculated the linear and nonlinear constitutive material relations in a moving medium. As expected from the transformation of electromagnetic fields, the movement of the medium leads to a mixing of the electric and magnetic response, where the obtained expressions for the linear constitutive equations were in agreement with the relations given in the literature. In addition, we showed that the movement of the medium also implies a mixing of the nonlinear material properties which describe cross-coupled nonlinear interactions between the E-and B-fields and the medium. This means, for example, that a Pockels medium at rest can display a Faraday effect if the material is moved relative to the observer, and vice versa.
In the final part of the paper we focused on time-independent, spatial coordinate transformations which are of particular interest for the design of nonlinear TO devices. As an illustrative example of such a device, we presented a twisted nonlinear field concentrator and calculated the second harmonic wave that is generated when the concentrator is illuminated by a strong laser field. In this respect, we demonstrated that sophisticated nonlinear phenomena in complex media can take a much simpler form if an appropriate coordinate transformation is applied.
The considerations have shown that the incorporation of nonlinear susceptibilities in the TO approach provides a promising computation method for calculating nonlinear effects in moving or inhomogeneous media and offers new opportunities for the design of novel optical devices with tailored nonlinear properties.

Appendix
Summary of the linear terms in equation (12) The linear terms in (12) (the free-space contribution and the linear response of the material) can be written as And finally, the last line of (12) follows after the re-definition 4πχ σκ µν (new) := 1 2 δ σκ µν − δ σκ νµ + 4πχ σκ µν .
Proof of equation (17) By using the identity the components of the electric displacement field are and, accordingly, the components of the magnetic fields are given by By comparing the expressions of D i and H i with we obtain the relations given in (17).

Number of independent components in
(exemplary for n = 2), there are only 6 independent combinations for each upper index pair α i β i , e.g. α i β i = 01, 02, 03, 12, 13, 23. Furthermore, in the sequence α 1 β 1 · · · α n β n all n pairs (α i β i ) can be permuted as a whole without changing the value of χ α 1 β 1 ···αnβn µν . Consequently, for fixed lower index pair µν, the number N of independent tensor components is equal to the number of combinations of n labels where each label can take 6 values (with repetition). This number is given by: For instance, N(1) = 6, N(2) = 21, N(2) = 56, etc. Finally, we can also vary the lower index pair µν for which also 6 independent combinations are possible due to the above symmetries. Hence, the total number of independent components of the tensor χ α 1 β 1 ···αnβn µν is given by 6N(n).
Intermediate transformations used in equation (45) The first intermediate transformation f is given by with the corresponding Jacobian matrix and determinant: The inverse intermediate transformation is given by: The corresponding Jacobian matrix and determinant are: