Gouy and spatial-curvature-induced phase shifts of light in two-dimensional curved space

Gouy phase is the axial phase anomaly of converging light waves discovered over one century ago, and is so far widely studied in various systems. In this work, we have theoretically calculated Gouy phase of light beams in both paraxial and nonparaxial regime on two-dimensional curved surface by generalizing angular spectrum method. We find that curvature of surface will also introduce an extra phase shift, which is named as spatial curvature-induced (SCI) phase. The behaviors of both phase shifts are illustrated on two typical surfaces of revolution, circular truncated cone and spherical surface. Gouy phase evolves slower on surface with greater spatial curvature on circular truncated cone, which is however opposite on spherical surface, while SCI phase evolves faster with curvature on both surfaces. On circular truncated cone, both phase shifts approach to a limit value along propagation, which does not happen on spherical surface due to the existence of singularity on the pole. An interpretation is presented to explain this peculiar phenomenon. Finally we also provide the analytical expression of paraxial Gaussian beam on general SORs. By comparing the result with the exact method we find the analytical expression is valid under the approximation that beam waist and scale of surface are beyond order of wavelength. We expect this work will enhance the comprehension about the behavior of electromagnetic wave in curved space, and further contribute to the study of general relativity phenomena in laboratory.


Introduction
It has long been acknowledged that when a converging, monochromatic beam propagates through its focus, it undergoes an additional phase shift of 180°, compared with the plane wave with same frequency. This phase shift is known as Gouy phase, which is discovered by Gouy in 1891 [1]. Basically Gouy phase is considered to originate from the transverse spatial confinement, as formulated by Feng and Winful [2], while various explanations have been reported from different perspectives [3][4][5][6][7]. As a fundamental physical phenomenon, this phase anomaly exists not only in electromagnetic beams but also in all kinds of waves which can be focused, such as matter wave [8][9][10], acoustic pulse [11,12], etc. On account of its wide existence and significant consequence, Gouy phase has been extensively investigated and applied in multifarious systems, such as Airy beams [13], surface plasmon polaritons [14], radially polarized beams [15], generation of few-cycle laser pulse [16], use in terahertz microscopy [17], effect on propagation of vector beam [18], etc. Among them the simplest and most well-known example is Gaussian beam in the background of flat space [19], where Gouy phase can be analytically expressed with respect to propagation distance zas j = ( ) z z arctan waves [25], etc. Among them, investigation of optics on curved surface is a burgeoning attempt [26]. By discarding one spatial dimension, researchers are able to fabricate the geometrical structure of two-dimensional space straightforwardly, while existence of massive body is not necessary. Experimentally it can be realized by covering a thin layer of waveguide on surface [27]. In the latest decade various concepts have been reconsidered and reported, such as solitons [28], evolution of speckle pattern [29], spatially accelerating wave packets following nongeodesic trajectories [30,31], topological phases in curved space photonic lattices [32], phase and group velocity of wave packets [33], Wolf effect of light spectrum [34,35], etc. Specially, in a pioneering work [26], Schrödinger equation for linear propagation on curved space with constant Gaussian curvature is derived, which is in essence the wave equation under paraxial approximation, and one of whose solutions is the fundamental notion in optics, Gaussian beam.
Paraxial approximation, corresponding to small divergence angle (or narrow angular spectrum) of beam, is commonly employed in optics studies. However, under some circumstances, when beam width is subwavelength or beam divergence is large, paraxial approximation is no longer valid. Therefore in this work we utilize angular spectrum approach, which is applicable for both paraxial and nonparaxial cases. This paper is organized as follows. In section 2, we calculate the output field by applying the generalized angular spectrum method. We obtain the expression of Gouy phase on curved surface, and meanwhile introduce an extra phase induced by curvature of surface. In section 3, we illustrate the evolutions of both phase shifts along propagation on two surfaces of revolution, circular truncated cone and spherical surface, and find distinction of their behaviors on different surfaces. A possible explanation is come up with in section 4, where the effect of curved space on light beam is physically revealed. In section 5, we provide analytical expressions of both phase shifts under paraxial approximation. In section 6, we present our concluding remarks.

Basic theory
In this work, we exemplarily investigate on a special class of curved surface-surface of revolution (SOR), whose rotational and translational symmetries in transverse direction remarkably simplify the problem. When appropriately choosing a set of orthogonal curvilinear coordinates, as shown in figure 1(a), the two-dimensional space is described by metric where r is the radius of revolution (ROR), and r 0 is the ROR at initial plane z=0 . Here the dependence of r on z is also the function of parametrization, which determines the way how surface is embedded in three-dimensional space. Wave equation of a ξ-polarized, time harmonic field on such curved surface is (for detailed derivation, see [26,36]) is the covariant Laplacian, g is determinant of metric tensor g g, ij is the element of inverse matrix of g. Here g=r(z)/r 0 , g 11 =1, and = ( ) g r r z is generally deemed as geometrical potential and plays a vital role in Hamiltonian of particles constrained on surface [37]. Concretely speaking, H and K are extrinsic and intrinsic curvature defined as average and product of main curvatures κ 1 , κ 2 , respectively. For general SORs, geometrical potential -H K 2 is expressed as which, in principle, cannot be neglected when the scale of surface is comparable to wavelength, or when ROR along z direction varies very quickly. By applying ansatz . When transverse spatial confinement is extreme, i.e. l ⪅ w 0 , the range of x k k may extend over 1, as in (c), and may further lead these modes to be evanescent.
and (c) are calculated by MAF method, and the imaginary parts are not exhibited as they are ignorable compared with real parts. induced by the covariant Laplacian, can be eliminated, whereas an extra effective potential Following [38], the input field Φ in (ξ) can be represented in the form of angular spectrum, which is essentially a superposition of infinite plane waves with different propagation directions (i.e. distinct transverse wave numbers), as follows where f(k ξ , ξ, z) represents plane-wave base, and A(k ξ ) is the corresponding superposition coefficient (or weight factor). Accordingly the output electric field is , with the subscript f being abbreviation of flat space. However, such plane waves are no longer solutions of wave equation on curved surface (i.e. equation (3)). Therefore, we should figure out the counterpart of plane wave on general SORs.
To this end, we follow the procedure in flat space and write , . By using the conventional method of variable separation, one can find that the transverse component c exp i , which is similar to that in flat space. Nevertheless, different from flat space where k ξ can be arbitrary continuous values, the geometrical structure of surface yields a periodical constraint that c , 2 0 should be fulfilled for arbitrary ξ, leading k ξ to be a series of discrete values m/r 0 , with m being an integer and indicating different discrete components. Consequently the integral in equations (4) and (5) should be substituted by summation. On the other hand, the longitudinal part is a solution of equation can be seen as the z component of effective wave number in curved space. Specially, when , this component becomes imaginary, and it correspondingly turns into evanescent. It is interesting to find that equation (6) is in the same form as a time-independent Schrödinger equation with effective potential varying over propagation distance z. Various approximation methods have been developed to solve equation (6), while here we adopt a combination of two representative approaches, Wentzel-Kramers-Brillouin (WKB) approximation and modified Airy function (MAF) method. WKB approximation is the most extensively used among all the approximation methods owing to its simplicity in derivation and reasonable accuracy. When approximated to the first order, Unfortunately, as WKB approximation is based on the assumption of G , inaccuracy and divergence occur in the vicinity of the so-called turning point . For components with considerable transverse wave number, or when parameters of surface satisfy some requirements (the requirements, however, vary from surface to surface), it is probable that turning point appears during propagation (i.e. > x z 0 k 0, ), and then WKB approximation breaks down. On this occasion we resort to MAF method which is fairly accurate both near and at turning point, and solution is given as (for detailed derivation of this method, see [39,40]) where Ai and Bi are Airy function of first and second kind, respectively, while C and Dare unknown constants which are going be determined as per boundary condition of specific surface, as will be shown below. The lower limit of integral in the expression of s x ( ) k z , is tactically chosen to be turning point so that this solution is exact when G ( ) z 2 is linear in z. Up to now we have generalized the concept of plane wave on arbitrary SORs. In view of the orthoganality of these plane-wave bases 1 , we ultimately come to angular spectrum representation of initial profile on general SOR And the output field, revised after equation (5), is ,e x p i .
Typical spectra of transverse wave number in source plane are illustrated in figures 1(b) and (c). The interval of transverse wave number between two modes is proportional to 1/r 0 , while width of spectrum increases with decrease of w 0 . So in practice, number of plane-wave components which contribute to the output field is determined by both r 0 and w 0 . For instance, in figure 1(b) where incident beam width w 0 is wide compared with wavelength λ, only few components should be counted because the mode interval is determined by r 0 . In contrast, in figure 1(c), when w 0 is narrower than λ, numerous components are involved in calculation, and compared with the width of spectrum envelop, the interval between two components is so small that summations in equations (10) and (12) can even be approximated to intergrals. By definition, Gouy phase is a finite-beam width-induced phase shift, and thus can be formulized as the difference between phase of output field and that of the (nondiffracted) plane wave propagating along longitudinal direction [6,13] is the ideal plane wave in 2D curved SOR. Besides, compared with plane wave propagating in flat space, another phase shift shows up apart from Gouy phase, which is arised by distinct physical origin. This phase shift is introduced by the curved surface itself, and thus is denominated as spatial-curvature-induced phase (SCI phase), and is defined as difference between phase of plane wave propagating along z direction in curved space and that in flat space, i.e.
From its expression, the SCI phase is merely subjected to longitudinal coordinate.

Circular truncated cone
We are going to apply the aforementioned theory to two typical SORs and investigate the properties of these two phase shifts. The first one is the circular truncated cone, which is fairly simple and can be parametrized by = + ( ) r z pz r 0 , where p describes the slope of surface. Turning point on such SOR is located at = -- , that is, when the mode is initially evanescent, and in the meantime the scale of surface is subwavelength. According to equations (7)-(9), one can obtain solutions as follow. For WKB approximation, in regions away from , e x pi , d . 1 These plane-wave bases satisfy relation that * 0 on both sides and integrating over ξ from p r 0 to π r 0 .
For MAF method, we have . Meanwhile in the region  ¥ z , WKB approximation is always valid, and the results obtained from these two methods should be unquestionably identical. By substituting Ai and Bi with their asymptotic forms [41], and then comparing it with equation (16), we finally have relation + = C D i 0. Without loss of generality, we take C=1, and = -D i. Figure 2 illustrates the behaviors of Gouy phase and SCI phase on such curved SORs, as well as the corresponding influence from initial ROR r 0 and beam waist w 0 . In figure 2(b), it is observed that along propagation, the absolute value of Gouy phase increases monotonically, and finally converges to a maximum. The maximum, as is shown in figures 2(a) and (b), increases with initial ROR r 0 , or namely, size of surface. The limiting case is when  ¥ r 0 , such SOR can be approximated to a flat space, with the corresponding Gouy phase being −π /4 for a 2D light field. Besides, when propagating on a certain SOR, the narrower beam will also lead to the greater maximal value of j | | G . The SCI phase, however, is independent on beam waist, as predicted in equation (15). From figures 2(d) and (e), SCI phase is always positive, which is different from Gouy phase. It increases with propagation distance and is saturated to a maximal value. More remarkably, this phase shift is more significant on the surface with smaller r 0 .

Spherical surface
Another classic SOR is spherical surface, which is the 2Dvisualization of a homogenous and isotropic universe with positive constant curvature [42]. Spherical surface is generally parametrized as , where R is the only constant to characterize such surface. Different from the circular truncated cone where the surface can extend to infinity, this model is closed and in effect two singularities exist at 'north pole' and 'south pole', as RORs there are infinitesimally small. Therefore without loss of generality, we restrict propagation distance z within [0, π R/2). In practice, the area quite close to singularity is abnegated since some unexpected phenomena caused by extremely small RORs and infinite potential, such as interference of beam with itself, may occur. The expression of output field can be obtained by WKB approximation and MAF method, following the same procedure performed on the circular truncated cone, except that here C=1, D=0. It is observed in figure 3(a) that compared with flat space, Gouy phase on spherical surface evolves faster; and the smaller the radius R, the faster Gouy phase evolves. Although due to constraint of limited propagation distance, it is unavailable to figure out maximal Gouy phase on spherical surfaces, tendency of the curves implies Gouy phase will not saturate to a fixed value as it behaves on circular truncated cone and in flat space. From figure 3(b), SCI phase increases during propagation, and smaller scale leads to greater phase shift, following the same rule as in last example. However, when approaching the singularity, SCI phase increases briskly and, theoretically speaking, will even approach infinity.

Discussion
In the following contents we are going to come up with a possible interpretation which may explain the above behaviors of Gouy phase and SCI phase. As formulated in [2], a finite transverse spatial extension of beam results in a spread of transverse momentum, which further reduces the axial component of wave vector (i.e. k z ). Gouy phase is originated from this deviation of k z from k. On arbitrary SORs, however, potential induced by surface, including geometrical potential -H K 2 and effective potential V eff (z), can be regarded as a part of wave number. For convenience, we can denote the summation of these two potentials as κ 2 (z), and consequently effective wave number . Now the curved surface can be equivalent to a medium with varying refractive index. Based on [2] and the Hamilton-Jacobi function, to which components of wave vector are subjected [26], the effective axial propagation constant is The inset shows the zoom-in of the area between two white dashed lines. Converging surface has a focusing effect on beam, so evolution of beam width is an interplay between surface and diffraction of beam itself. In this case beam continues being focused during propagation. In this figure, we choose WKB approximation in calculation. order to roughly estimate the validation of the paraxial approximation, we can define the parameter 0 as a criterion. In general, when s=1, the paraxial approximation is always valid. However, if κ 2 (z) decreases with z as shown in figure 2(c), the minimal value of k eff is given by . Therefore the criterion of the paraxial approximation here can be estimated simply by <  ( ) s kw 1 1 0 . When κ 2 (z) increases monotonically with z, this condition may become Figure 4 takes circular truncated cone as example and illustrates to what extent these analytical expressions accord with exact solutions. It can be observed that when beam waist w 0 and initial ROR r 0 are subwavelength, there is conspicuous difference between these two methods. However, when beam waist and scale of surface increase to or are beyond order of several wavelengths, these analytical expressions can already serve as an eligible approximation of exact results. In figures 4(e) and (f), they show evolution of Gouy phase and SCI phase on curved surfaces with different curvature as well as in the flat space. Under paraxial approximation, when r 0 increases, the potential induced by surface decreases, that is, effect of curvature is weaker. In the limit of zero curvature, evolution of Gouy phase approaches to the situation of flat space, while SCI phase approaches to zero.

Conclusion
In conclusion, compared with flat space, curvatures (including extrinsic and intrinsic curvature, as well as the effective potential which stems from transformation in covariant Laplacian and is essentially related to intrinsic curvature [27]) of two-dimensional curved space introduce an effective potential, enabling curved surface to be equivalently treated as an imhomogenous media, and thus induce an extra phase shift. Another influence of curved SOR is that the proper length varies at different position because of the variation of RORs. This may explain the acceleration and deceleration of Gouy phase on curved SORs. These properties may offer a feasible and unprecedented avenue for controlling phase shifts of light fields in curved space. Moreover, The generalized angular spectrum method we develop is applicable in both paraxial and nonparaxial regime, which may also result in a potential practicability in exploring numerous optical phenomena in curved space.