A Unified Approach to Nonlinear Transformation Materials

The advances in geometric approaches to optical devices due to transformation optics has led to the development of cloaks, concentrators, and other devices. It has also been shown that transformation optics can be used to gravitational fields from general relativity. However, the technique is currently constrained to linear devices, as a consistent approach to nonlinearity (including both the case of a nonlinear background medium and a nonlinear transformation) remains an open question. Here we show that nonlinearity can be incorporated into transformation optics in a consistent way. We use this to illustrate a number of novel effects, including cloaking an optical soliton, modeling nonlinear solutions to Einstein’s field equations, controlling transport in a Debye solid, and developing a set of constitutive to relations for relativistic cloaks in arbitrary nonlinear backgrounds.

Transformation optics [1][2][3][4][5][6][7][8][9] , which uses geometric coordinate transformations derive the materials requirements of arbitrary devices, is a powerful technique. Essentially, for any geometry there corresponds a material with identical transport. With the correct geometry, it is possible to construct optical cloaks [10][11][12] and concentrators 13 as well as analogues of these devices for other waves [14][15][16][17][18] and even for diffusion [19][20][21][22][23][24][25][26] . While many interpretations and formalisms of transformation optics exist, such as Jacobian transformations 4 , scattering matrices [27][28][29][30][31] , and conformal mappings 3 , one of the most theoretically powerful interpretations comes from the metric formalism 32 . All of these approaches agree that materials define an effective geometry, however the metric formalism is important since it allows us to further interpret the geometry. In particular, certain geometries correspond to solutions to Einstein's field equations, which relate geometric curvature to gravitational forces. Materials that mimic these geometries, or artificial relativistic media, constitute a subset of transformation optics materials (dark blue circle, Fig. 1) that can effectively model relativistic effects 32 , such as black holes 33,34 and gravitational lensing 35 or create novel devices such as the space-time cloak (which hides events instead of objects) 36 .
One limitation of transformation optics, however, is the necessity of using linear materials (materials whose properties do not change with electric field, pressure, temperature, etc.). At present, the transformations that have been derived have exclusively been applied to linear media. That is, the focus has been upon media equivalent to an isotropic, homogeneous, linear background medium embedded in curvilinear coordinates. However, there is no necessity to maintain the constraint of linearity. In thermal transformations, researchers have already considered the case of temperature dependent transformations (which we shall generalize as "nonlinear transformations"), and shown how they are equivalent to a thermally nonlinear material embedded within a linear background 37,38 . However, considerations of background nonlinearity have thus far been absent. Moreover, nonlinear transformations lack the intuitive physical interpretation of linear transformation materials, where transport follows stationary geodesics. This intuition is useful when developing devices where geodesics bend and shift with the applied field.
In this paper, we shall present a unified theory of nonlinear transformation optics. We will consider both the case of a nonlinear background medium (bottom half of Fig. 1) and nonlinear transformations (right half of Fig. 1) in arbitrary combination. We shall begin by generalizing transformation optics theory to incorporate these nonlinearities, then consider examples illustrating this formalism from each of the new, nonlinear quadrants of Fig. 1. Examples will be selected for their practical significance, physical intuition, and clarity.

Results
Nonlinear Transformation Formalism. To begin, in linear transformation optics, the constitutive relation is 32 where γ ij is the metric of the reference frame (unprimed coordinates) and we have used Einstein summation notation for curvilinear coordinates (indices repeated as both subscript and superscript (covariant and contravariant) are summed, Latin indices are only over spatial dimensions, Greek indices are over space and time (0 th ) dimensions). Note that this formulation explicitly requires a stationary medium and observer. A fully covariant description is given in [39][40][41][42] and is easily modified to incorporate the formalism developed here, but is less widely used than the stationary case.
Scalar field dependence. Since in most cases the variable of interest in a transformation materials problem is a scalar (temperature, pressure, electric potential) or can be approximated as such (electric field, fluid velocity, etc. for specific geometries), it is helpful to begin the consideration of field dependence by considering this scalar variable explicitly. By definition, this scalar variable is not changed by transformation media techniques except for the implicit change of spatial variable (e.g. T(r) → T(r′)). Moreover, scalars remain scalar even when raised to an arbitrary power, so rank of the tensors considered in equation 1 are not changed by the field dependence (even when they are Taylor expanded as a function of field strength). As the tensor rank governs the rules for the transformation, we conclude that equation 1 can easily be generalized to a nonlinear transformation of a nonlinear background by the relation where we have assumed that the nonlinearity is solely a function of implicitly time and coordinate-dependent and position (henceforth an implicit depedence, with dependence upon other fields (e.g. magnetic field intensity, temperature, pressure, etc.) following trivially from this equation.

Figure 1.
Representation of our transformation optics framework. Background material () and coordinate transform (g) can be linear or nonlinear with respect to applied fields, making four mutually exclusive cases. Within this parameter space, certain combinations satisfy transformation optics requirements (light blue ellipse). A subset of these also satisfy Einstein's field equations (dark blue ellipse). When nonlinearity is included, effects from other fields, e.g. nonlinear optics (magenta circle) can become incorporated into transformation optics. Transformation optics is typically fabricated using metamaterials (orange circle), although other implementations are possible (if often trivial, e.g. lenses) and applications of metamaterials outside of transformation optics also exist. Examples of transformation optics devices from each quadrant are labelled, with the nonlinear examples being explored in the text (except the switchabe cloak, discussed in ref. 37 ). Note that the functional forms of E ( )  and g(E) are arbitrary. Assuming the coordinate transformation x i → x i′ leaves Maxwell's equations (or the corresponding equation of motion for other fields) unchanged, except for a change of variables (i.e. [ ( ), ( , ( )), ( ( )), ] 0 0 for operator  that defines E), then the introduction of nonlinearity preserves transformation optics techniques, as the coordinates only enter the nonlinearity through the field. For example, this means that the Kerr nonlinearity (discussed in detail in the following section) should not be perfectly transformed according to the general rules for coordinate transformations. In particular, because that nonlinearity depends upon the electric field j which is changes form under coordinate transformations, Maxwell's equations are changed if the transformation is applied blindly 43,44 . If this effect is neglected, however, Maxwell's equations actually obey the correspondence of equation 4 and ′ E x E x ( ) ( ), as desired.
Tensor field dependence. This scalar field nonlinearity is a useful form but is not the most general nonlinearity. In general, the field of interest will be an arbitrary rank tensor. However, if we assume that the functional dependence of g,  are accurately represented by Taylor expansions in terms of field strength, then we need only consider the special case of a vector field depedence. Higher ranked tensors will result in terms which are effectively reducible to tensor products of vectors (e.g. ). Because the rank of the tensor determines the properties under coordinate transformation, we thus see that the rules for the vector field dependent transformation cover the rules for any other tensor rank (assuming that the Taylor expansion holds). In particular, a material with anisotropic dependence upon the electric field (i.e. both strength and orientation) can be represented by the Taylor expansion (2 ) is the nonlinear (electric) susceptibility (magnetic susceptibility can be defined similarly for μ(H)). Since each factor of E i in the anisotropy corresponds to an additional scalar product, where indices ij are handled by the transformation of equation 3 which multiplies all factors of . Note that this formalism can implicitly handle anisotropy in indices ij, as such linear anisotropy can be introduced as a transformation from an isotropic background, giving a composite transformation.
Nonlinear Background -Linear Transform. In particular, if the nonlinearity takes the form of a Kerr nonlinearity    (P is polarization and χ is susceptibility, which we assume to be isotropic), Maxwell's equations remain unchanged under the cloaking transformation, Thus, if we can find a solution to Maxwell's equations in Euclidean space with a Kerr nonlinearity, we can find a solution to Maxwell's equations with a Kerr cloak permittivity (lower left in Fig. 1) by writing the Euclidean solution in primed (i.e. cloak) coordinates (note that we have taken advantage of the transformation E(r) = E(r′) and assumed that only one polarization of electric field is present 43 ). (Including the full anisotropy of χ (3) , equation 9 is modified χ χ , more details are given in 43,44 ) The Kerr nonlinearity is a special case of nonlinear optics with an exactly solvable system for special values of intensity E 2 corresponding to optical soliton modes. For concreteness, we select the first spatial soliton 45 , is the pulse width, k is the wave-vector, ω frequency, (3) 0 2 0 2 , and n 0 is the linear index of refraction. Note that the polarization and direction of soliton propagation in this case are arbitrary, but because the equation 9 is isotropic with respect to these choices the system does not need to be modified for solitons of different polarization or propagation direction. (Higher orde solitons, however, will possess different combinations of frequencies and so will require some modification in practice for proper realization of the cloak.) The analytic solution is plotted in Fig. 2. Note that the cloaking is exact in the analytic case, despite the nonlinear background. While this solution is exact, it is limited to planar system (i.e. 1 + 1D) such as slab waveguides, since the Kerr spatial soliton is unstable in higher dimensions 46 . However, this approach is easily generalized to other systems with optical spatial solitons, including saturable Kerr nonlinearity, photoreactive materials, etc. Additionally, as discussed below, the spectral breadth of the soliton implies that the cloak cannot work perfectly for all the frequencies present (but can work perfectly for some and approximately for others, as the soliton spectrum is relatively narrow).
Applications. Moreover, while the ability to cloak a soliton is somewhat artificial, its success implies that an observer could not use a background Kerr nonlinearity to detect the presence of a nonlinear cloak, whereas a linear cloak would disperse the soliton and thereby render itself detectable. This could have applications in controlling waves in other nonlinear media with similar nonlinearities, such as protecting against rogue waves 47 (although this is easier to accomplish in optics than surface waves, as the nonlinearity is more easily engineered there). In addition, nonlinearity is used to locate defects in acoustic/ultrasonic nondestructive testing [48][49][50] , including via acoustic solitons 51 , and to detect waves in sonar systems 52 . Applying this technique to an acoustic system, then, would prevent an observer from using the nonlinearity of a material to detect a cloak hidden within it. Since nonlinearities exist in most real materials (e.g. polycrystalline solids, porous media like clay or soil, and seawater and other complex fluids), the incorporation of nonlinearity to cloaking design is integral to rigorously prevent detection of buried objects or evade sonar. Similarly, a concentrator (equation 30) could be used to increase the density of a field in a region, thereby promoting soliton formation or promoting the interaction of multiple solitons.
Realization. Given that materials with a Kerr nonlinearity exist, it is likely that the simplest implementation of a soliton cloak is merely to incorporate such a material into existing metamaterial cloak designs. However, such cloak designs are based upon resonant effects (e.g. superluminal phase velocity at a single frequency) and therefore cannot be directly applied to multiple disparate frequencies. The lowest order spatial soliton, though, actually preserves the spectal properties of a plane wave fairly well. Compared to a plane wave in vacuum, there are only two modifications: the wave-vector along zˆ is shifted from k to k − γ (i.e. still a single frequency) and additional components are introduced along k y according to the Fourier transform of sech, which is sech(πy 0 k y /2). This is an approximately Gaussian spread centered at zero and with standard deviation 2/πy 0 . While this implies that the soliton only possesses frequencies in the neighborhood of the corresponding plane wave frequency, the Gaussian peak is an issue for perfect cloak operation. It is not, however, an impediment to imperfect cloaking. If, instead of mapping a finite region to a point, the scattering cross-section of a domain is reduced by a factor Δr, then the resulting cloak is subject to the constraint where Δω is the bandwidth of the cloak's operating frequency 53 . If the operating frequency is greater than soliton's spectral width, then the approximate cloak still works and is in principle experimentally feasible.
Linear Background -Nonlinear Transform. While we have seen that transformation optics is robust to background nonlinearity, that case is easier to understand. The dynamics there are identical to nonlinear optics in Euclidean space, with the added linear transformation merely distorting the geodesics in fixed directions. When the transformation is nonlinear, then the geodesics can change with changing intensity. This makes, say, the combination of a nonlinear transform and the Kerr effect far harder to calculate. Instead, we shall now consider only a nonlinear transform and fix the background to be linear (upper right in Fig. 1). We can apply physical intuition to the nonlinear transform by taking inspiration from the study of effective gravitational fields via linear transformation optics 33,34 , where variations in the permittivity mimic the gravitational field produced by a mass distribution. In that case, the metric used must satisfy Einstein's field equations where R μν is the Ricci curvature tensor , , , c is the speed of light, and T μν is the stress-energy tensor. That is, a matter distribution is used to define T μν , which then defines g μν via Eq. 12, thereby defining the equivalent , μ via Eq. 1. However, relativity also predicts that energy and mass are equivalent (as in the famous = mc 2  ). As such, energy distributions can also define a stress-energy tensor and thereby produce a gravitational field 54 .
If the only source of energy is the electromagnetic field, then solutions to Eq. 12 are referred to as electrovacuum solutions. A material satisfying Eq. 3 with a metric transform obeying Eq. 12, then, will have a nonlinearity equivalent to the gravitational field produced by the electromagnetic field.
For a purely electromagnetic source, T μν is ) . For reference, contrast this with observations of a local, inertial observer, who will define For both inertial and non-inertial observers with for a purely electromagnetic source, = μ μ T 0 so R = 0 and our equations simplify. Additionally, we shall use the intuitive interpretation of a local observer to denote T 00 with U(E), simplifying our notation.
To be more specific, we consider a plane wave solution (note that other solutions exist for more realistic fields than a monochromatic plane wave, e.g. 55 , but are not as easily expressed in terms of the effective material parameters needed for transformation materials). In particular, our electromagnetic field is selected such that a local observer will detect a plane wave, which leads to 2 where f is defined below in terms of the metric, u = ω(t − x/c). If our background is linear, then = U E E ( ) 0 2  (note that this is not the averaged energy, it retains space and time dependence), S cU E x ( ) → =ˆ, and σ δδ x . We can then assume a metric of the form g μν ≡ diag[−1, 1, f(ct − x), f(ct − x)] in Minkowski coordinates and use Eq. 12 to get defining h = f 2 . Using the identity 2 cos ( ) 1 cos(2 ) 2 φ φ = + , the stability condition h → 1 as |E| → 0, and the rotational symmetry (implying −|E| should give the same solution as |E|), gives where MathieuC is the Mathieu cosine function (which, because the first term is negative, behaves closer to cosh than cosine) and E c G / 4  is the natural electric field scale. Note that G only occurs in E 0 , and so an effective gravitational effect can be tuned by changing E 0 . In Fig. 3a, we plot f, where we've exploited the periodicity of Eq. 17 to create a periodic continuation of f (using the unmodified form results in an exponential growth of ). We now consider the composite transform T GR,C = T C T GR , to create a cloaked region within this artificial relativistic medium. Using f, and Eqs 3 and 8 we calculate the components , , xx xy yy    and plot them in Fig. 3b-d. Notably, we do not plot E for this setup, as it is indistinguishable from the solution to a purely linear cloak. This is expected, given that the form of E was assumed in solving for g, but we can also show that Maxwell's equations reduce to  Nonlinear Background -Nonlinear Transform. While the Mathieu cosine form of the nonlinear transform is helpful for illustrating the physical relevance of a nonlinear transform, an alternative formulation is preferable for developing materials prescriptions. In particular, it is preferable in nonlinear optics to know the dependence of the susceptibility as a Taylor series in E to some finite order. In considering this problem, we shall once allow  E ( ) to have an arbitrary nonlinearity, as that is most useful for design (bottom right of Fig. 1). Eq. 17 remains unchanged, save for a modification of U(E) to reflect the new value of E ( )  , but the solutions can no longer be expressed in terms of analytic functions. Instead, we employ the Liouville-Neumann series technique to solve for h(u), where u = ω(t − x/c). That is, we consider a series expansion h E h u ( ) We truncate our solution f = h 2 at 4 th order in E, as terms of that order and below are most relevant to nonlinear optics. However, truncation means that our solution takes the form f(E = |E| sin(u), u), as some terms have a more depend upon u m (i.e. secular terms from nonlinear resonance) that cannot be factored without higher order terms (these likely correspond to the cosh dependence in the Mathieu solution). This solution f gives the vacuum nonlinearity in a flat space-time. We now apply the cloaking transformation Eq. 8 in cylindrical coordinate to this metric and use Eq. 3 to get in transformed cylindrical coordinates (μ ij /μ(E) defined identically). To verify that this leaves Maxwell's equations unchanged, it suffices to observe that this can also be written (the nonlinear resonance should be modified to functions of x′ as they derive from E(r)), so  clk ( ) will reproduce Eq. 8, transforming Maxwell's equations from a set of operators L r r as desired for a nonlinear cloaking transformation. Notice that the first term in Eq. 22 is the standard linear cloak (recall that f(E = 0) = 1 and the second is purely due to the vacuum nonlinearity. We can thus define Notice that, although E 0 was originally defined in terms of the constant G, it is the only place that such constant enters into the transformed material equation. Thus, we are free to redefine E 0 as any effective scale for the electric field strength, rather than the scale prescribed by Eq. 12. That is, we can use transformation optics to model a nonlinear gravitational field with arbitrary strength E 0 (G eff ).

Realization. Equation 23
includes nonlinear effects up to fourth order, which is considerably higher than most applications for nonlinear optics. Although effects like four wave mixing explicitly rely upon these fourth order nonlinearities (hence their inclusion in equation 23), it is far more common to retain only the lowest order non-trivial nonlinearity. Assuming our background lacks an  (1) term, which is trivially decoupled from the gravitational nonlinearity, the first non-trivial term is the secondorder term: However, the kl dependence of χ (3) is left implicit here and the u 2 dependence gives a non-stationary nonlinearity. This issues can be corrected to lowest order by making the approximation u ≈ sin u to get |E| 2 u 2 − E 2 ≈ |E| 2 − 2E 2 . If we interpret the direction of our input field as fixed, this gives ( ) On the other hand, recalling that our original choice of the direction of E was arbitrary within an isotropic background gives us Either form of χ gives the simplest model of the gravitational nonlinearity, which can be seen as a more complex form of Kerr nonlinearity. Applications. Given the complexity present in even the simplest possible realizations of this gravitational nonlinearity, it is unlikely that these nonlinear gravitational transformations will find experimental applications in the near future. However, these results hold an important theoretical application -namely the elucidation of the correspondence between nonlinear optical media and electrovacuum gravitational fields. This implies that difficult problems in general relativity can be mapped to equivalent problems in nonlinear optics, where approximations as to the effect of nonlinearity are more robust. Simultaneously, results in nonlinear optics can now be interpreted as effectively gravitational, providing an interesting intuitive picture of photon-photon interactions not normally considered in the standard treatment.
Transformation Media Extension. Before considering our final example, it is worth stepping back and considering how these nonlinear transformation optics techniques could be extended to other forms of transformation media. Acoustics is by far the easiest generalization, as there are straightforward mappings from transformation optics to transformation acoustics 17 . Heat transport and diffusion are more difficult, however. While the introduction of field dependence to the already established thermal transformation 19 holds -i.e. that −is valid, the diffusion equation is not Lorentz invariant and therefore is not a valid equation for the relativistic interpretation. Thus, while transformation materials is applicable to nonlinear heat transport, it cannot be interpreted in terms of effective gravitational fields. However, because transformation diffusion is defined for an isotropic background κ 0 , with all anisotropy arising from the transformation, a further interpretation is plausible. Both the background nonlinearity and isotropic nonlinear transform control the speed of diffusion at a given temperature, whereas the anisotropic aspect controls the preferential direction of diffusion as a function of temperature. So, for our final example we consider heat transport within a Debye solid (κ ∝ (T/T 0 ) 3 , c p ∝ (T/T 0 ) 3 , ρ = ρ 0 / (1 + αT) ≈ ρ 0 , where T 0 is the Debye temperature, and α is thermal expansivity (O(10 −5 /K for a solid)). As the nonlinear transform in this case is an arbitrary g(T) that does not satisfy Eq. 12, we consider a "phase transition" transform where λ = ρc p ,κ rr or κ θθ , λ L(H) are the low (high) temperature transformed parameters, T tr is the transition temperature, and T Δ is the range of the intermediate zone.
In particular, we want a cloak for high temperatures (Eq. 8, T > T tr ) and a concentrator for low temperatures. For a concentrator's transformation a central r < R 2 is shrunken down to r < R 1 where R 1 < R 2 . A second region, R 2 < r < R 3 is then stretched out to accommodate the shrunken central region. Contrast with equation 8, where a single point is expanded to a finite diameter, leaving a hole in the center and only requiring one transformed region a < r < b. Typically 20 , the concentrator uses a radially linear mapping   , which is helpful since many discretized thermal metamaterial implementations of the thermal cloak/concentrator rely upon this symmetry being preserved 19,20 . It is thus conceivable to consider building such thermal metamaterials necessary for the cloak-concentrator out of a binary composite where one material undergoes a phase transition at a specific temperature. For concreteness, we define the overall form of the cloak-concentrator's material parameters (including both spatial and thermal dependance in the annular region) Note that for T T tr  equation 34 converges to equation 31, for  T T tr equation 34 converges to equation 32, and T ≈ T tr is equal to the average of equations 31 and 32 but switches quickly to either the high or low temperature limit when T T T tr | − | > ∼ Δ . This behavior is confirmed by COMSOL simulations, where the far field temperature distribution in the nonlinear background is unchanged by the presence of the transformation materials (Fig. 4a) and that the low (Fig. 4b) or high (Fig. 4c) temperature cases work as a thermal concentrator or cloak of the Debye solid. More interesting is the intermediate case when T(x = 0) ≡ T tr ≈ 8.4T 0 , i.e. when transition temperature isotherm bisects the annulus. In this regime the device acts like a cloak for < x 0  and a concentrator for > ∼ x 0, Fig. 4d. A careful examination of x ≈ 0 reveals that neither the device is neither a cloak or a concentrator close to the transition temperature, with isotherms shifting to reflect this compromise regime. In particular, the thermal gradient within the central domain neither vanishes (as in the cloak) or a constant (as in the concentrator, note the isotherms in Fig. 4b), but curved towards the low temperature phase. Moreover, within cloak-concentrator there is a clear asymmetry to the slope of the isotherms for x > 0 vs x < 0 in the high and low temperature phases but not in the compromise region (the slope of some isotherms even oscillates near x = 0). However, this intermediate region could be made extremely narrow in practice, as it corresonds to a coexistance of the high and low temperature phases of the thermal metamaterial and conventional materials generally have unambiguous state for first order phase transition.
Applications. With this functionality we have created a device which excludes heat flux at high temperatures and attracts it at low temperatures, thereby acting as an effective thermal regulator. In addition, there is ability of the cloak (equation 31) and concetrator (equation 32) designs to function correctly when embedded within temperature dependent materials is a non-trivial improvement of the functionality of thermal transformation materials, as all realistic materials will possess some temperature dependence to their parameters. It is also worth noting that even physical parameters unrelated to thermal transport will generally possess some form of temperature dependence, and therefore that the incorporation of designs for transformation materials with non-trivial temperature dependent backgrounds or transformations could greatly improve the versatility of transformation Scientific RepoRts | (2018) 8:4436 | DOI:10.1038/s41598-018-22215-x materials. Similar arguments can be made for materials with dependence upon pressure, electric, magnetic, and electromagnetic fields, albeit for narrower ranges of materials.
Realization. Since most materials possess some degree of temperature dependence to their thermal properties, the realization of the cloak-concentrator seems fairly straightforward and is perhaps the most easily achieved of all the examples discussed in this work. To work through how it might be realized, however, we shall make several approximations. In particular, as our goal is the simplest experimental realization, we shall consider the design of the bilayer cloak 21 as our starting point. This design, which used a pair of annular regions of homogeneous, isotropic material to create a thermal cloak, explicitly focused upon the steady state. This implied that the heat capacity could be completely neglected, as engineering the conducitivity was all that was required. More subtly, the bilayer cloak also assumed that the superposition holds in its derivation, and therefore only works rigorously for materials without temperature dependence. As such, we shall neglect the temperature dependence of the background and assume that the materials in our cloak-concentrator are approximately independent of temperature except for near their respective (and possibly distinct) transition temperatures, where the operation of the device breaks down. Making these assumptions turns the problem into a two-fold one: first a bilayer concentrator must be designed to match the implementation of the bilayer cloak. And second there must be a pair of materials which can, after undergoing phase transitions, meet the thermal conductivity requirements in each layer. To address the first part, consider the logic of the bilayer cloak -to minimize the thermal flux in the center the inner ring's conducitivity was set as close to zero as possible and the outer ring's conducitivity was then adjusted to prevent any scattering at the outer edge of the cloak. For the concentrator a similar requirement holds, the outer layer must prevent scattering, but the inner layer's conductivity must be adjusted to maximize the heat flux in the interior (which is filled with the same material as the background). Numerically optimizing the analytic solution to Laplace's equation with these constraints gives properties as laid out in Table 1, which enhances the heat flux in the concentrator regime by 8 percent. Note that the temperature dependence of the conductivities in each layer are opposed to each other, with the inner layer's conducitivity falling with temperature and the outer layer's conductivity rising with rising temperature. This is a non-trivial temperature dependence, but not outside of the realm of experimental possiblity -notably thermal diodes have been implemented using pairs of materials with similar thermal properties. As such, the thermal cloak-concentrator is likely the ideal candidate for the experimental verification of nonlinear thermal transformations, particularly if the constraint of κ i → 0 for the cloak is relaxed (i.e. an approximate cloak design is used.

Conclusions
In summary, we have developed a formalism for understanding transformation materials in nonlinear media and undergoing nonlinear transformations, and shown how this formalism can be applied to soliton transport (with applications in hiding objects from nonlinear detection protocols or shielding rogue waves), effective gravitational fields (with applications in modeling relativistic effects in nonlinear optical materials), and thermal management (with applications in controlling heat flow in realistic, thermally nonlinear materials or incorporating temperature dependent parameters in arbitrary transformation media devices). It is therefore possible to use nonlinear transformation media to model a wider variety of transport phenomena than have been previously considered. Furthermore, the constitutive relations that we have derived in equations 3, 23, 27 and 34 can be used for a wider variety of transformations than just the cloaking transformation (e.g. concentrators, rotators, camoflague, etc.).

Additional Applications.
Given the incorporation of frequency dependent material parameters in the metamaterial realizations of transformation optics 11 and the development of nonlinear metamaterials 56 , the incorporation of nonlinear optical transformations presents the opportunity for a wider variety of functional materials in this framework. Wave mixing, for example, could be used to generate additional field components at a frequency several times the incident wave's. If a metamaterial has different resonant effects at these frequencies, it could (say), screen the incident wave while concentrating the nonlinear contribution. This could be accomplished with metamaterials whose parameters are quite similar to equation 34, except that the dependence upon (T − T tr )/T Δ would go to a dependence upon (ω − ω tr )/ω Δ . Or it could change the relative field concentration within different regions (e.g. via a concentrator or its inverse, a dilutor), thereby controlling the effective strength of the local nonlinearity (increasing the effective nonlinearity in high fields regions and decreasing it in low field regions), which could alternatively aide in the observation of novel nonlinear effects or prevent the occurence of unwanted nonlinearities. This approach is quite similar to that applied by certain classes of nonlinear metamaterials, which rely upon concentrating the electromagnetic field on regions of strong nonlinearity 57 or improving the efficiency of nonlinear conversion (e.g. via phase matching) 58 . Soliton formation 59 could also be promoted, say by counteracting an excessive dispersion or nonlinearity in the external medium (an alternative approach to that adopted by 60 ).

General Realizations.
In principle, any of these applications could be implemented with the utilization of nonlinear metamaterials, which have so far developed independently of advances in transformation materials, being more focused upon the creation of arbitrary strength nonlinearities. Nonlinear metamaterials can be produced by adding nonlinear elements to a linear metamaterials design 61 or embedding such designs in nonlinear backgrounds 62 , or introducing intrinsic structural nonlinearities to linear materials (e.g. magnetoelastic structural reconfiguration) 63 . Such nonlinear metamaterials have already been used for elementary soliton steering applications 60 . Note that all of these techniques are reviewed in more detail in Lapine et al. 56 and the references therein. The combination of frequency mixing and resonance is not necessarily a fatal flaw to these designs either. While most of the devices considered here (cloaks, concentrators, etc.) rely upon a tailored metamaterial resonance for their implementation, the potential existence of waves at multiple, distant frequencies can be accounted for. In particular, as the metamaterial resonance is tied to the length scale of the structuring, the incorporation of multiple length scale would allow for additional resonances. Thus, a cloak designed for a frequency doubling application should be designed with two resonant lengt scales (the doubled frequency corresponding to half the length). While that complicates the fabrication, it should be achievable so long as the specification application is well-defined in advance. Generally, nonlinear transformations can increase the versatility of frequency-dependent phenomena in metamaterials and transformation materials.  Table 1. Properties of the temperature independent cloak-concentrator. Cloak parameters taken from ref. 21 , with concentrator conductivity chosen to optimize heat flux given the geometry.