Light-driven mass density wave dynamics in optical fibers

We have recently developed the mass-polariton (MP) theory of light to describe the light propagation in transparent bulk materials [Phys. Rev. A 95, 063850 (2017)]. The MP theory is general as it is based on the covariance principle and the fundamental conservation laws of nature. Therefore, it can be applied also to nonhomogeneous and dispersive materials. In this work, we apply the MP theory of light to describe propagation of light in step-index circular waveguides. We study the eigenmodes of the electric and magnetic fields in a waveguide and use these modes to calculate the optical force density, which is used in the optoelastic continuum dynamics (OCD) to simulate the dynamics of medium atoms in the waveguide. We show that the total momentum and angular momentum in the waveguide are carried by a coupled state of the field and the medium. In particular, we focus in the dynamics of atoms, which has not been covered in previous theories that consider only field dynamics in waveguides. We also study the elastic waves generated in the waveguide during the relaxation following from atomic displacements in the waveguide.


I. INTRODUCTION
The angular momentum of light has recently gained much interest both in free-space and fiber-based optical communications [1][2][3][4][5][6][7][8][9] as it has been suggested that the orbital angular momentum of light could provide an additional degree of freedom for data multiplexing [1,2].The angular momentum of light has been investigated in the quantum domain [10,11], in various photonic materials [12][13][14][15], and also in the near-field regime [16][17][18].Also encoding/decoding using superpositions of spatial modes in an optical fiber has been recently demonstrated [19].Despite the long history of the studies of optical waveguides [20,21], the detailed theoretical description of the angular momentum of light in waveguides has been missing.In particular, the polariton property of light as a coupled state of the field and the medium and its relation to the linear and angular momenta of light have been considered for dielectric materials only recently [22][23][24].
Combining the electrodynamics of continuous media, elasticity theory, and Newtonian mechanics, we have shown that a light pulse propagating in a medium drives forwards an atomic mass density wave (MDW), which carries a major part of the total linear and angular momenta of light in many common dielectrics [22][23][24].In the classical optoelastic continuum dynamics (OCD), formulated by us, the MDW is an inevitable consequence of the optical force density of the electromagnetic field.When the optical force density and the elastic forces are included in the Newtonian equation of motion, one can solve numerically the space-time dynamics of the medium.By considering the numerically calculated space-time dynamics of the atoms in the medium, one can see that field intensity maxima are accompanied by atomic density maxima, which propagate with the same velocity as the field intensity maxima.Note, in particular, that these MDWs should not be confused with elastic waves, the dynamics of which is governed by elastic forces and which propagate with the velocity of sound.Instead, the MDWs should sooner be considered as shock waves, which are generated by the optical force field.
In the single-photon picture, the coupling of the electromagnetic field to the atomic MDW gives rise to masspolariton (MP) quasiparticles, which are covariant coupled states of the field and matter having a nonzero rest mass [22,23].In this MP quasiparticle model, the MDW follows directly from the Lorentz transformation of the energy and momentum of light.Hence, the MP is fundamentally different from the conventional excitonpolariton and the phonon-polariton quasiparticles, in which a photon is in resonance or in close resonance with an internal excited state of the medium.Within the framework of the MP theory of light, the correct way of thinking the exciton-polariton and the phonon-polariton quasiparticles is to think them as coupled states of the MP and an exciton or an optical phonon, respectively.Note, however, that in the present work we limit our investigations to highly transparent materials in which the coupling of the MPs to the excitons or optical phonons is very weak, and therefore, we do not elaborate this coupling further in this work.
The MP quasiparticles have been shown to carry a total momentum of the Minkowski form p MP = n p ω/c, where n p is the phase refractive index of the medium, is the reduced Planck constant, ω is the angular frequency, and c is the speed of light in vacuum [23].The total MP momentum is split between the electromagnetic field and the MDW so that the share of the field corresponds to the Abraham momentum p field = ω/(n g c), where n g is the group refractive index.The MDW carries the difference of the Minkowski and Abraham momenta.Therefore, the MP theory of light provides a win-win resolution to the centennial Abraham-Minkowski controversy of pho-ton momentum in a medium [25][26][27][28][29].
In this work, we apply the MP theory to simulate the propagation of the MDW of selected light pulses in a step-index circular waveguide.It is well-known that, in waveguides, part of the electromagnetic field propagates as an evanescent field in the cladding layer surrounding the core.In the MP theory of light, this affects the atomic MDW that propagates partly in the core and partly in the cladding layer.The total magnitude of the the MDW is, however, fully determined by the covariant state of light in a medium [22,23].If the cladding layer is vacuum, the entire MDW is propagating in the core material and its magnitude should be reduced when compared to the MDW in the case of a homogeneous dielectric.This reduction of the MDW should also be related to the waveguide dispersion.Here we show that the numerical OCD simulations accurately verify these points and the full correspondence with the MP quasiparticle model.This also shows that the optical force density used in the OCD model in homogeneous materials also accurately applies to waveguides.
Note that there exist some previous works that address the question of consistent definition of the angular momentum of light in circular waveguides [30][31][32].These works do not, however, consider the essential role of the atomic MDW as a carrier of the major part of the total angular momentum of light in waveguides.Therefore, these works have inherent limitations like (1) not being able to predict the transfer of atomic density with light, (2) the inevitable dissipation of field energy resulting from the transfer of atomic density, and (3) these theories considering only fields cannot account for losses and dephasing effects that result from the coupling of the MDWs with other fundamental excitations of the medium.In contrast, these aspects are either directly obtained from the numerical OCD simulations or can be calculated using the properties of the coupled fieldmedium MP state.
This work is organized as follows: Section II reviews the theoretical foundations of the OCD model.Section III presents the solution of the electric and magnetic fields of a light pulse in a step-index circular waveguide geometry.The OCD simulations of the propagation of selected light pulses in a circular waveguide are presented in Sec.IV, where we also describe the elastic relaxation dynamics of the nonequilibrium atomic displacements that the MDW leaves in the waveguide.In Sec.V, we compare the OCD results with the results obtained by using the MP quasiparticle model.Finally, conclusions are drawn in Sec.VI.

A. Newton's equation of motion
Light is optoelastically coupled to the medium by Newton's equation of motion.This is one of the key elements of the OCD model.Since the force density of the optical field can give only small velocities for the atoms of a continuous medium, the atomic velocities are certainly nonrelativistic.Therefore, Newton's equation of motion for the mass density of the medium ρ a (r, t) can be written as where r a (r, t) is the instantaneous position-and timedependent atomic displacement field of the medium, f opt (r, t) is the optical force density experienced by the medium atoms, and f el (r, t) is the elastic force density, which arises between atoms that are displaced from their initial equilibrium positions by the optical force density.
In previous literature, several forms of the optical force density have been extensively discussed [33,34].We have recently shown, for a low-loss dielectric, there is only one form of optical force density that is fully consistent with the MP quasiparticle model and the underlying principles of the special theory of relativity [22].This force density is presented for nondispersive media as [22,33] This optical force density leads to the transfer of an inevitable part of the total momentum of light by the MDW [22,23].In this work, we show that the optical force density in Eq. ( 2) is also applicable in the description of optical forces of confined modes, such as light propagating in optical fibers or waveguides.
In the OCD simulations, we use the well-known elastic force density that follows from Hooke's law [35].Using the material displacement field r a (r, t), this elastic force density can be given for an isotropic linear elastic medium as [36] (3) where λ L and µ L are the Lamé elastic constants of the medium.The Lamé elastic constants could also be replaced by any two independent elastic moduli, such as the bulk and shear moduli whose relation to the Lamé constants is described in Ref. [37].For anisotropic cubic crystals, such as silicon, the elastic force density has a slightly more complicated form, which is given, e.g., in Refs.[23,35].
When we apply the OCD model to simulate the propagation of light pulses in optical waveguides in Sec.IV, we use a perturbative approach in which we neglect the extremely small damping of the field that results from the transfer of field energy to the kinetic and elastic energies of the medium atoms by the optical force density in Eq. (2).In the case of light propagation in a homogeneous medium, the accuracy of this approximation is estimated in Ref. [22], but the conclusions are also valid for waveguides provided that optical losses are small, which is typically a good approximation as light can travel in waveguides over long distances.
B. Energies, momenta, and angular momenta of the field and the mass density wave The total electromagnetic energy of the light pulse is obtained as an integral of the instantaneous energy density of the field and the total mass energy of the MDW is obtained as an integral of the MDW mass density just as in the case of a homogeneous medium [22].These relations and the expression of the total MP energy of a light pulse are given by As shown in Ref. [22], the kinetic energy of the MDW is extremely small and neglected in E MDW , which only describes the mass energy of the MDW.Again, in the same way as in a homogeneous medium, the total momentum of the coupled MP state of the field and matter and the momentum shares of the electromagnetic field and the atomic MDW are given by [22] In the MP theory of light, we also obtain the total angular momentum of the MP, denoted by J MP , and its shares J field and J MDW carried by the electromagnetic field and the atomic MDW.These angular momenta are given by [24]

III. ELECTRIC AND MAGNETIC FIELDS OF A LIGHT PULSE IN A CIRCULAR WAVEGUIDE
We aim at using the electric and magnetic fields of a light pulse in a waveguide to simulate the propagation of the MDW driven by the light pulse.Therefore, we next review the solution of the electric and magnetic fields of a light pulse in a step-index circular waveguide geometry illustrated in Fig. 1.The waveguide core radius is R. Inside the core, the refractive index is n 1 while the refractive index in the cladding layer is n 2 .

FIG. 1.
Illustration of the step-index circular waveguide geometry.The waveguide core with radius R has a refractive index n1 while the refractive index of the cladding layer is n2.The total length of the waveguide is L.

A. Fields in a step-index circular waveguide
The electric and magnetic field modes are given as exact solutions of the Maxwell's equations for the stepindex circular waveguide geometry as presented, e.g., in Ref. [38].These solutions are called Bessel modes.The necessary condition for the existence of confined modes is n 1 > n 2 .The cladding layer is assumed to be thick enough so that the field of confined modes is virtually zero in the outermost boundary of the cladding layer and we only need to account for the core and the cladding when computing the spatial distribution of the electromagnetic field inside the waveguide.
Here, we assume a light pulse with a nonzero spectral width as described by the Gaussian function u(k 0 ) = e −[(k 0 −k0)/∆k0] 2 /2 /( √ 2π∆k 0 ), where k 0 = ω 0 /c is the wavenumber in vacuum for central frequency ω 0 , and ∆k 0 is the standard deviation of the wavenumber in vacuum.In the cylindrical coordinate system with coordinates r, φ, and z and basis vectors r, φ, and ẑ, the instantaneous electric and magnetic fields of the Bessel mode pulses are given, in the waveguide core (r < R), as and, in the cladding layer (r > R), as Here J l (x) are the Bessel functions of the first kind, K l (x) are the modified Bessel functions of the second kind, , and ω(k z ) = ck z /n p,eff is the waveguide dispersion relation with the effective phase index denoted by n p,eff and the wavenumbers of the core and cladding given by k 1 = n 1 k 0 and k 2 = n 2 k 0 , respectively.The derivatives of the Bessel functions are denoted by J l (x) and K l (x).The other primed quantities in Eqs. ( 7) and ( 8) correspond to the values obtained for k 0 instead of k 0 .The materials are here assumed to be nondispersive.The unknown amplitude factors A, B, C, and D can be related to each other by applying the boundary conditions at the interface between the core and the cladding as described below or, e.g., in Ref. [38].
In this work, we apply the fields above to describe the propagation of light pulses with a relatively small spectral width.In this limit, the light pulse is approximately monochromatic and Eqs. ( 7)-( 10) can be approximated further in the same way as presented for Laguerre-Gaussian pulses in Ref. [24].

B. Characteristic equation and the relations between the field amplitudes
The unknown amplitude factors A, B, C, and D can be related to each other by the boundary conditions.Requiring that E φ , E z , H φ , and H z are continuous at the interface between the core and the cladding at r = R leads to four equations, given by Equations ( 13) and ( 14) lead to the characteristic equation for k z from which the coefficients A and B have been cancelled.This characteristic equation is transcendental and it reads For a given l and a given frequency ω, there will be only a finite number of solutions k z that satisfy Eq. ( 15).Once a specific solution is known, one can use Eqs.( 11)-( 14) to solve the ratios of the amplitude factors A, B, C, and D. The absolute magnitudes of the amplitude factors become determined independent of a free phase factor by normalizing the fields so that the total electromagnetic energy of the light pulse is E field .

C. Waveguide dispersion
The effective phase refractive index n p,eff for the field propagating in the optical fiber is given as the ratio of k z , solved from the characteristic equation in Eq. ( 15), and k 0 as [38] In terms of the total MP momentum in Eq. ( 5) and the field energy in Eq. ( 4), the effective phase refractive index can also be written as n p,eff = c|P MP |/E field .From the characteristic equation, one also obtains the effective group refractive index n g,eff , given by There also exists an alternative well-known expression for the group refractive index as a ratio of the integrals of the energy density and the Poynting vector.Using Eqs. ( 4) and ( 5), we then obtain n g,eff = E field /(c|P field |).

TE and TM modes
First, we consider the special cases of the T E and T M modes.For the T E modes, we have E z = 0, from which it follows that A = C = 0.In this case, the left hand sides of Eqs. ( 13) and ( 14) become zero.Equation ( 14) is then satisfied if l = 0 and Eq. ( 13) leads to the characteristic equation for the T E modes given by where we have used the Bessel function relations J 0 (x) = −J 1 (x) and K 0 (x) = −K 1 (x).The remaining Eq. ( 11) is trivially satisfied and Eq. ( 12) relates the nonzero field amplitudes B and D to each other.
Respectively, for the T M modes we have H z = 0, from which it follows that B = D = 0.In this case, the right hand sides of Eqs. ( 13) and ( 14) become zero.Equation ( 13) is then satisfied if l = 0 and Eq. ( 14) leads to the characteristic equation for the T M modes given by where we have again used the Bessel function relations J 0 (x) = −J 1 (x) and K 0 (x) = −K 1 (x).The remaining Eq. ( 12) is trivially satisfied and Eq. ( 11) relates the nonzero field amplitudes A and C to each other.
The T E 0,m and T M 0,m modes are numbered so that the lowest mode number m corresponds to the solution of the characteristic equation for which the value of hR is the smallest.The cutoff condition for the existence of the T E 0,m and T M 0,m modes is given by [38] where x 0,m , m ≥ 1, is the mth zero of the Bessel function J 0 (x).The first three zeros of the Bessel function J 0 (x) are x 0,1 = 2.405, x 0,2 = 5.520, and x 0,3 = 8.654.

EH and HE modes
Second, we consider the general solution of the characteristic equation in Eq. (15).Solving the characteristic equation for J l (hR)/(hRJ l (hR)) gives Here, the plus sign corresponds to the the EH modes and the minus sign corresponds to the HE modes.The cutoff conditions for the existence of the the EH 1,m and HE 1,m modes are given by [38] where x 1,m , m ≥ 1, is the mth zero of the Bessel function J 1 (x) and x 1,0 = 0.The first three zeros of the Bessel function J 1 (x) are x 1,1 = 3.832, x 1,2 = 7.016, and x 1,3 = 10.173.The cutoff values of (R/λ 0 ) 1,m for the HE 1,m modes are smaller than for the EH 1,m modes, and especially, the mode HE 1,1 does not have a cutoff.The cutoff conditions for the existence of the the EH l,m and HE l,m modes with l ≥ 2 are given by [38] where x l,m , m ≥ 1, is the mth zero of the Bessel function J l (x) and z l,m is the mth root of equation

E. Angular momentum and relation to Laguerre-Gaussian modes
The Bessel modes studied in this work are eigenmodes of step-index circular waveguides as discussed above.The step-index profile is used in most singlemode fibers and in some multi-mode fibers.In counter distinction, Laguerre-Gaussian modes are eigenmodes of graded-index circular waveguides with a parabolic refractive index profile.The parabolic profile is the most common profile in graded-index multi-mode fibers as it minimizes modal dispersion by continual refocusing of the rays in the fiber core.The Bessel and Laguerre-Gaussian beams can both carry angular momentum.This has been experimentally demonstrated by rotating microparticles in optical tweezers [39][40][41][42][43][44][45][46].

IV. SIMULATIONS
Next, we apply the OCD model to illustrate the node structure of the atomic MDW and the actual atomic displacements due to a light pulse in a circular waveguide made of fused silica.The light pulse of Eqs. ( 7)-( 10) is assumed to have a total electromagnetic energy of U 0 = 1 µJ and a vacuum wavelength of λ 0 = 1550 nm.The corresponding photon number is N ph = E field / ω 0 = 7.803 × 10 12 .In our examples, the waveguide core radius is R = λ 0 /2.The relative spectral width of the pulse used in the simulations is ∆ω/ω 0 = ∆k 0 /k 0 = 0.05, which corresponds to the temporal full width at half maximum (FWHM) of ∆t FWHM = 27 fs.The FWHM is fixed to this small value, which is close to the engineering feasibility limit, to make the node structure of the MDW visible in the same scale with the longitudinal Gaussian envelope of the pulse.In our simulations, we use temporal discretization of h t = 2π/(40ω 0 ) and spatial discretization of h x = λ/40, where λ = λ 0 /n p,eff is the effective wavelength of the studied mode in the waveguide.This discretization is sufficiently dense compared to the scale of the harmonic cycle of the optical field.The flowchart of the OCD simulation is described in Appendix C of Ref. [22].
The simulations are performed for a waveguide whose core is assumed to be made of fused silica, which is the standard material used in optical fibers.The cladding layer is assumed to be vacuum meaning that there is no cladding.The corresponding fields are a special case of Eqs. ( 7)- (10) obtained by setting n 2 = 1.The phase and group refractive indices of fused silica are given by n p = 1.4440, and n g = 1.4626 for λ 0 = 1550 nm [47].As the material dispersion is not very large, in this work, we neglect the material dispersion and use the phase refractive index only.However, the material dispersion could be accounted for as explained in Ref. [23].The waveguide dispersion is here the dominating form of dispersion and it is described by the effective phase and group refractive indices of the waveguide geometry following from Eqs. ( 16) and (17).For the T E 0,1 mode studied below, the effective phase refractive index of our waveguide geometry is n p,eff = 1.1330 and the effective group refractive index is n g,eff = 1.5856.Respectively, for the HE 1,1 mode that is studied in the second example, the effective phase refractive index is n p,eff = 1.3033 and the effective group refractive index is n g,eff = 1.5355.The mass density of fused silica used in the OCD simulations is ρ 0 = 2200 kg/m 3 [48], and the Lamé elastic constants are λ L = 15.4GPa and µ L = 31.3GPa [49].

A. Simulations of the atomic mass density wave
First, we study the atomic MDW of a temporally Gaussian T E 0,1 light pulse for which l = 0. See Media 1 for a video file of the simulation.Figure 2(a) shows the longitudinal atomic velocities of the simulated MDW as a function of position in the plane y = 0 µm when the T E 0,1 pulse center is propagating at the position z = 0 µm.In vacuum, outside the waveguide core at x < −λ 0 /2 or at x > λ 0 /2, the atomic velocities are zero as the atomic MDW does not exist.Inside the waveguide core, the atomic velocities in the MDW follow the field intensity.The node structure of the field is then also well seen in the atomic MDW.The atomic velocities are zero at the center of the waveguide at the position x = 0 µm following the T E 0,1 mode profile.When the Gaussian envelope of the pulse drops to zero on the left of z = −5 µm and on the right of z = 5 µm, also the atomic velocities of the MDW become zero.
Figure 2(b) shows the longitudinal atomic displacements of the simulated MDW as a function of position along z axis in the plane y = 0 µm when the T E 0,1 pulse center is propagating at the position z = 0 µm corresponding to the atomic velocities in Fig. 2(a).One can see that the atomic displacements are largest at lateral positions of maximum field intensity.At positions z > 5 µm, the atomic displacement is zero as the light pulse has not yet reached these positions, and therefore, these atoms have not yet experienced the optical force density.
In vacuum, outside the waveguide core, the atomic displacement in Fig. 2(b) is zero as the atomic velocities in Fig. 2(a).After the light pulse at z < −5 µm, the atomic displacement has obtained a value that is in the femtosecond time scale practically constant in time, but which varies in the lateral direction.The longitudinal atomic displacements in the waveguide cross section just after the light pulse has gone is presented in Fig. 2(c).These atomic displacements correspond to the atomic displacements after the light pulse in Fig. 2(b).The atomic displacement averaged over the waveguide cross section is given by ∆r MDW = δM/(ρ 0 πR 2 ) = 2.135 × 10 −15 m.The atomic displacements that the light pulse leaves behind are relaxed to equilibrium only in the much larger time scale of elastic waves.The relaxation dynamics of the waveguide in the longer time scale is studied in more detail below.
The OCD simulation results in Fig. 2 are fully consistent with the MP quasiparticle model of Refs.[22,23] as the total transferred mass and momentum following from the atomic velocities and displacements equal the quasiparticle model results within the numerical accuracy of the simulation.For the total transferred mass, the MP quasiparticle model gives an expression δM = (n p,eff n g,eff − 1)N ph ω 0 /c 2 that is obtained from the homogeneous field result presented in Ref. [23] by simply replacing the phase and group refractive indices of the material with the conventional effective values for the studied mode in the waveguide geometry.The numerical OCD model gives the total transferred mass for our T E 0,1 mode light pulse as δM = [ρ a (r, t) − ρ 0 ]d 3 r = 8.862 × 10 −24 kg, which corresponds to the value obtained by using the MP quasiparticle model formula.For the momentum of the MDW, the OCD simulations give p MDW = | ρ a (r, t)v a (r, t)d 3 r| = 1.676 × 10 −15 kgm/s, which corresponds to the MP quasiparticle model value obtained by using p MDW = (n g,eff − 1/n g,eff )N ph ω 0 /c.
Another interesting case is the HE 1,1 mode light pulse for which l = 1.This pulse does not have a cutoff, and in contrast to the T E 0,1 pulse above, it carries angular momentum of per photon.See Media 2 for a video file of the simulation.Figure 3(a) depicts the atomic velocities of the simulated MDW as a function of position in the plane y = 0 µm when the HE 1,1 pulse center is propagating at the position z = 0 µm.Inside the waveguide core, the atomic velocities in the MDW again follow the field intensity.As the HE 1,1 mode is not a transverse mode, the pulse does not have as clear peaks and throughs as the T E 0,1 pulse above.Also, in contrast to the T E 0,1 pulse, the atomic velocities obtain their maximum values at the center of the waveguide at the position x = 0 µm in accordance with the HE 1,1 mode profile.On the left and right of the light pulse, the atomic velocities in the MDW are again zero as expected.
Figure 3(b) depicts the longitudinal atomic displacements of the simulated MDW as a function of position along z axis in the plane y = 0 µm when the HE 1,1 pulse center is propagating at the position z = 0 µm.These atomic displacements correspond to the atomic velocities in Fig. 3(a).One can again see that the longitudinal atomic displacements are largest at lateral positions of maximum field intensity.Before the light pulse at positions z > 5 µm, the atomic displacement is again zero and, after the light pulse at z < −5 µm, the atomic displacement has obtained a value that is temporally approximatively constant in the femtosecond time scale.The longitudinal atomic displacements of the HE 1,1 pulse just after the light pulse has gone is presented in the waveguide cross section in Fig. 3(c).These atomic displacements correspond to the atomic displacements after the light pulse in Fig. 3(b).The atomic displacement averaged over the waveguide cross section is given by ∆r MDW = δM/(ρ 0 πR 2 ) = 2.683 × 10 −15 m.
The simulation results are again in full correspondence with the MP quasiparticle model of Refs.[22] and [23].For the total transferred mass of our HE 1,1 pulse, we obtain δM = 1.114 × 10 −23 kg and, for the momentum of the MDW, we obtain a value p MDW = 2.175 × 10 −15 kgm/s.These values are slightly larger than the values in the case of the T E 0,1 pulse due to the smaller evanescent field that propagates in vacuum outside the fiber core and that does not directly contribute to the atomic MDW.In the MP quasiparticle model, this difference between the quantities is accounted for the effective phase and group refractive indices of the waveguide corresponding to the somewhat smaller waveguide dispersion in the HE 1,1 case.Figure 3(d) presents the transverse components of the atomic velocities in the MDW at a particular instance of time in the middle of the HE 1,1 light pulse.The transverse components of the atomic velocities follow from the nonzero transverse component of the optical force density and the Poynting vector, which does not exist in the case of the T E 0,1 pulse above.As described in Ref. [24], the transverse components of the atomic velocities are related to the transfer of a substantial part of the total angular momentum of light in a medium.The angular momentum of the electromagnetic field of the HE 1,1 pulse obtained using the OCD model is J field = 0.586N ph ẑ and the angular momentum of the MDW is J MDW = 0.414N ph ẑ.The sum of these angular momentum contributions is the total angular momentum of the coupled state of the field and matter, given by J MP = N ph ẑ, which is essentially an integer multiple of .This result verifies the entanglement of the field and medium components of the total angular momentum of a light pulse in the case of a waveguide geometry.Previously, this entanglement has been studied in the case of a homogeneous medium in Ref. [24].
The cross section of the total transverse atomic displacement field that the HE 1,1 pulse leaves behind, and which follows from the transverse atomic velocities in Fig. 3(d), is depicted in Fig. 3(e).Together with the longitudinal atomic displacements in Fig. 3, these transverse atomic displacements are again relaxed to equilibrium in the time scale of elastic waves.

B. Simulations of the elastic relaxation dynamics
Next, we simulate the beginning of the relaxation of the nonequilibrium atomic displacements of the waveguide after the T E 0,1 light pulse of Fig. 2 is gone.The initial atomic displacement for this relaxation process is that presented in Fig. 2(c).We assume that the fiber is infinitely long, which is realized in the simulations by setting periodic boundary conditions in the longitudinal direction just after the light pulse is gone.The damping of the elastic waves in a closed system is not accounted for by the elastic force density in Eq. ( 3).Therefore, the elastic waves that we obtain do not attenuate.Hence, we cannot simulate how the equilibrium of atomic positions is obtained at the end of the relaxation.Instead, we focus on the first nanosecond when the damping of the elastic waves can be neglected.
Figure 4 presents the relaxation dynamics of the waveguide after the T E 0,1 light pulse of Fig. 2 is gone.It shows the longitudinal atomic displacement in the waveguide cross section at three instances of time in the picosecond time scale.See Media 3 for a video file of the simulation.Figure 4(a) shows the longitudinal atomic displacement at t = 110 ps, Fig. 4(b) at t = 200 ps, and Fig. 4(c) at t = 260 ps.One can see that the minima and maxima of the atomic displacement switch their positions in the course of time.In Figs.4(a) and 4(c), the maximum atomic displacement is located at the origin while, in Fig. 4(b), at the origin, the atomic displacement obtains its minimum value.This can be understood as elastic shear waves that propagate from the circular fiber boundary to the center and back.The atomic displacement minima and maxima alternate, but their time-evolution is not fully periodic since the initial state of the relaxation process is not an eigenstate of the elastic oscillation of the waveguide.
Figure 4(d) shows the kinetic and strain energies of the elastic relaxation waves per unit length of the fiber as a function of time.At t = 0, the elastic energy is pure strain energy as the light pulse has only displaced atoms from their initial equilibrium positions but does not have left any net kinetic energy to them.At t > 0, the kinetic and strain energies alternate as a function of time so that their sum remains constant due to the conservation of energy.The highest peaks and throughs in the kinetic and strain energies at t = 413 ps and t = 826 ps are seen to be located close to times at which the shear waves with velocity v ⊥ = µ L /ρ 0 = 3770 m/s have propagated distances that are multiples of the waveguide diameter.
The value of the total elastic energy can be used to estimate the effective imaginary part of the refractive index n i , which corresponds to the dissipation of the field energy when the field propagates in the waveguide.The total elastic energy per unit length is 1.1 × 10 −18 J/m as can be read from Fig. Dividing the dissipated by the field energy U 0 of the pulse, gives the attenuation coefficient 1.1 × 10 −12 1/m.When this is compared with the conventional formula the attenuation coefficient, given by α = 2n i k 0 , we obtain n i = 1.4 × 10 −19 .This is vastly smaller than the imaginary part of the refractive index due to other physical nonidealities in highly transparent optical waveguides.
Since the velocity of light vastly larger than velocity of the process practically only after our short light pulse is gone.In this approximation, during the relaxation process, the total external force on the medium atoms is zero.Therefore, in accordance with Newton's first law, the center of mass of the waveguide remains at rest since the light pulse has only displaced atoms from their initial equilibrium positions but does not have left any net kinetic energy to them as discussed above.
The same kind of relaxation waves, as obtained for the longitudinal atomic displacements in Fig. 4, are also obtained for the transverse atomic displacements in the relaxation of the atomic displacements due to optical modes carrying angular momentum.These waves are also elastic shear waves in nature.Compressional waves are not obtained as the light pulse does not leave any mass density differences behind in the middle of the waveguide.As shown in Ref. [22], compressional waves are found in the case that a light pulse crosses an interface between materials of different refractive indices.Then, the optical surface force displaces atoms at the surface and these displacements are relaxed in terms of both compressional and shear waves.In the case of optical fibers, compressional waves can originate from the processes at the ends of the fiber or at any positions where there are changes in the refractive index.

V. COMPARISON OF THE OCD AND MP QUASIPARTICLE MODELS
The correspondence between the OCD and MP quasiparticle model results is summarized in Table I.The re-TABLE I.The transferred mass, the total momentum, the field's share of the momentum, and the MDW's share of the momentum calculated by using the MP and OCD models.The MP column shows the per photon value obtained by dividing the total quantities of the MP model with the photon number of the pulse N ph = E field / ω0.
sults for the energy and momentum of the MP and its field and medium parts are identical to the results for a homogeneous medium in Ref. [22] with an exception that we must now use the effective phase and group refractive indices of the waveguide in the MP quasiparticle model.Consequently, we obtain the ratio n p,eff n g,eff − 1 for the MDW's and the electromagnetic field's shares of energies and momenta.This result is directly related to the Lorentz transformation and the covariance principle of the special theory of relativity in the same way as described in the case of a homogeneous medium in Ref. [22].
The waveguide dispersion also influences the separation of the total angular momentum of light between the field and matter as described by the parameter ζ = |J field |/|J MP | used in Table I.Analytic calculation using Eq. ( 6) and the fields in Eqs. ( 7)- (10) shows that the angular momenta of the field and the MDW are given by  In the case of a homogeneous medium with a refractive index n, we have ζ = 1/n 2 as shown in Ref. [24].In contrast, in the waveguide geometry studied in this work, we have not found a simple expression for ζ in terms of the effective phase and group refractive indices of the waveguide.One can argue that this is not surprising as the covariance condition of the special theory of relativity used in the MP quasiparticle model relates the energies and momenta of the field and the MDW to each other [22] but no such a simple quasiparticle model relation is known to exist for the separation of angular momentum.In the case of angular momentum, only the total angular momentum of which is a multiple of , is well known to independent of the inertial reference frame.

VI. CONCLUSIONS
In conclusion, we have used the MP theory of light related OCD model to simulate the atomic MDW of propagating light pulses in a step-index circular waveguide geometry.We found unambiguous correspondence between the MP quasiparticle and OCD models in the description of light propagation in waveguides.Since part of the total energy and momentum of light propagates in vacuum outside the waveguide core, the magnitude of the MDW is This reduction place a way that is fully accounted for by using the well-known effective phase and group refractive indices of the waveguide in place of the corresponding indices of a homogeneous medium in the MP quasiparticle theory for dispersive media formulated by us recently.
In this work, we have also studied the elastic relaxation dynamics of the atomic displacements that the light pulse has left in the waveguide.We have found that these displacements are relaxed by elastic shear waves.The OCD model allows detailed modeling of the optoelastic dynamics of any material geometries under the influence of the optical force field, e.g., one can simulate the MDW in typical waveguides used in engineering applications.The OCD simulations can also be used for the optimization of possible measurement setups planned for the experimental verification of the atomic MDW associated with light.
After the present work was completed, the preprint in Ref. [50] including results complementary to our work, came to our attention.This preprint does not, however, separate the total angular momentum of light into the field and the atomic MDW contributions as we do in the present work.Therefore, the coupled dynamical description of the electromagnetic field and the atoms in the waveguide is not covered in the mentioned work.

5 FIG. 2 .
FIG.2.Simulation of the mass transfer due to a temporally Gaussian T E0,1 mode light pulse in a silica fiber (see alsoMedia  1).(a) Simulated longitudinal atomic velocity of the MDW in the fiber as a function of position in the plane y = 0 µm.(b) Longitudinal atomic displacement in the fiber as a function of position in the plane y = 0 µm.(c) Longitudinal atomic displacement in the fiber cross section z = 0 µm just after the pulse has gone.The pulse energy is E field = 1 µJ, the wavelength is λ0 = 1550 nm, and the duration is ∆tFWHM = 27 fs.The radius of the circular fiber core is R = λ0/2 and it is surrounded by vacuum.The dashed line shows the boundaries of the fiber.

1 FIG. 3 .
FIG. 3. Simulation of the mass transfer due to a temporally Gaussian HE1,1 mode light pulse in a silica fiber (see also Media 2).(a) Longitudinal atomic velocity of the MDW in the fiber as a function of position in the plane y = 0 µm.(b) Longitudinal atomic displacement in the fiber as a function of position in the plane y = 0 µm.(c) Longitudinal atomic displacement in the fiber cross section z = 0 µm just after the pulse has gone.(d) Instantaneous transverse atomic velocity of the MDW in the fiber as a function of position in the fiber cross section at the position of the pulse center at z = 0 µm.(e) Transverse atomic displacement in the fiber as a function of position in the fiber cross section z = 0 µm just after the light pulse is gone.The pulse energy is E field = 1 µJ, the wavelength is λ0 = 1550 nm, and the duration is ∆tFWHM = 27 fs.The radius of the circular fiber core is R = λ0/2 and it is surrounded by vacuum.The dashed line shows the boundaries of the fiber.

1 FIG. 4 .
FIG.4.Simulation of the elastic relaxation dynamics of the silica waveguide after the T E0,1 mode light pulse (see alsoMedia  3).The pulse energy is E field = 1 µJ, the wavelength is λ0 = 1550 nm, and ∆tFWHM = 27 fs.Longitudinal atomic velocity of the MDW at z = 0 µm cross section is shown at (a) t = 110 ps, (b) t = 200 ps, and (c) t = 260 ps after the pulse has gone.The radius of the circular fiber core is R = λ0/2 and it is surrounded by vacuum.The dashed line shows the fiber boundary.(d) The alternation of the elastic (kinetic and strain) energies per unit length of the fiber as a function of time.The elastic waves are gradually damped and thermalized by absorption and scattering not included in our simulations.Due to the periodic boundary conditions used in the simulation of the elastic relaxation, the interface effects at the ends of the fiber are not present.