Microscopic analysis of the energy, momentum and spin distributions in a surface plasmon-polariton wave

We analyze the electromagnetic field near a plane interface between a conductive and a dielectric media, under conditions supporting surface plasmon-polariton (SPP) propagation. The conductive medium is described by the hydrodynamic electron-gas model that enables a consistent analysis of the field-induced variations of the electron density and velocity at the interface and its nearest vicinity. The distributions of electromagnetic dynamical characteristics: energy, energy flow, spin and momentum are calculated analytically and illustrated numerically, employing silver-vacuum interface as an example. A set of the"field"and material contributions to the energy, spin and momentum are explicitly identified and classified with respect to their physical origins and properties, and the orbital (canonical) and spin (Belinfante) momentum constituents are separately examined. In this context, a procedure for the spin-orbital momentum decomposition in the presence of free charges is proposed and substantiated. The microscopic results agree with the known phenomenological data but additionally show specific nanoscale structures in the near-interface behavior of the SPP energy and momentum, which can be deliberately created, controlled and used in nanotechnology applications.


Introduction
Localized light fields attract significant attention due to their impressive abilities within optical nano-probing, precise optical manipulation and optical data processing. Especially, the surface plasmon-polariton (SPP) waves emerging near the interface between dielectric and conductive media are intensively investigated in connection with nano-and micro-optical means for light guiding, switching and controlling, which is crucial for further microminiaturizing of optical information devices and systems [1][2][3][4][5][6]. Many stimulating and promising applications are related with unique dynamical properties of the SPP waves such as transverse spin and momentum [7][8][9][10][11], special spin-momentum locking [11,12], nonreciprocity and unidirectional propagation [11][12][13][14].
The exclusive properties as well as the very existence of the SPP waves are directly associated with the dispersion (frequency dependence of the main electric and magnetic parameters) of both contacting media. Together with the strong spatial inhomogeneity of the SPP-supporting systems, this poses challenges in the theoretical description of the SPP field, especially its dynamical characteristics (DCs): energy, energy flow, momentum, angular momentum and their derivatives. The main problems relating the proper DC description in such "complex" media were recently resolved in a series of works [17][18][19][20][21][22][23]. This approach is, in fact, a generalization of the known Brillouin's method for the energy description in lossless dispersive media, which is based on the special dispersive corrections of their permittivities and permeabilities [24,25]. In application to the field momentum and angular momentum, it is closely related to the canonical (spin-orbital) momentum decomposition [26][27][28]: the dispersion-caused corrections are applied not to the usual Poynting momentum density but to its orbital (canonical) and spin (Belinfante) parts separately.
When being applied to the SPP field, this methodology offers a physically meaningful and consistent description of its momentum and angular momentum with their dependence on the media's properties and the radiation frequency [20][21][22][23]. Importantly, this phenomenological approach is justified by the microscopic analysis [20][21][22] based on the simplest lossless Drude model [2,5] for the electron gas in the metal. However, the Drude model entails some contradictory conclusions relating the SPP field spatial behavior in the near-interface region: the discontinuity of the transverse electric field and non-zero electron velocity at the metal surface lead to non-physical "singular" values of the SPP-induced electron density and to non-realistic data for the plasmonic field enhancement [29].
On the other hand, the microscopic justification of the phenomenological results for the spin and orbital SPP-momentum constituents raises a new problem in the electromagnetic DC theory. Actually, the existing results on the spin-orbital momentum decomposition essentially employ the dual symmetry [30] expressed by the "electric-magnetic democracy" [27] of the Maxwell equations without free charges and currents. In the microscopic electron-gas analysis, the presence of charges and currents destroys this symmetry, making the form and the very existence of a consistent algorithm for the momentum canonical decomposition unclear.
As we will see, both difficulties mentioned in the above paragraphs require a more detailed investigation of the DC spatial distributions with involvement of the electrons' quantum properties, and in the simplest way this can be performed on the ground of the hydrodynamic "jellium" plasma model [29,[31][32][33]. The hydrodynamic model has widely been used in the SPP theory; in particular, some steps were also undertaken in connection with the SPP momentum decomposition [22] but the systematic DC analysis was not accomplished. The solution of this problem is one of the present paper's aims.
In further Sections we present a comprehensive hydrodynamic microscopic description of the SPP field with the main attention paid to its DCs. The first objective is to justify the phenomenological results, including those for the spin-orbital (canonical) momentum decomposition; actually, this task was partly realized [20,22] but now we seek a more reliable record of the microscopic hydrodynamic corrections including those that are usually neglected in the phenomenological picture but are important for the physical consistency. In this process, we perform a second objective: a consecutive identification and classification of the meaningful "blocks" which compose the SPP momentum. This classification includes the "field" and material contributions, differing by their origin and physical nature, as well as their orbital (canonical) and spin parts. And finally, the obtained analytic results are used for the study of the DC properties, especially their spatial details associated with the hydrodynamic model and localized in the nearinterface region.
The presentation starts in Section 2 with reproducing the known solutions for the field and material characteristics within the hydrodynamic model framework. This part is intended for completeness and convenience; besides, the known results are presented in the forms suitable for subsequent comments and references. In Section 3 these solutions are applied for deriving the refined expressions of the field DCs. Especially, subsection 3.2 illustrates the calculation of the field-induced material momentum following the general scheme of Ref. [16] and finally presents the SPP momentum as a sum of the "field" and material contributions of different origins; subsection 3.3 describes the general algorithm for the spin-orbital momentum decomposition and details of its application to separate "field" and material contributions. The peculiar features of the SPP angular momentum (which is of the spin nature and has an "extraordinary" transverse direction [7][8][9][10]) are analyzed in subsection 3.4. In Section 3, the consideration is concentrated on the DC spatial distributions whereas the radiation frequency (the SPP wavelength) is kept constant; in contrast, Section 4 discusses the phenomena owing to the variable radiation frequency.
The meaning and peculiar features of the new hydrodynamic-based results are disclosed via the systematic juxtaposition with the usual phenomenological expressions for the SPP field and its DCs. To this end, the important known data, summarized and presented in the forms suitable for comparative discussions, are collected in the Supplemental Document. General significance of the main results, their possible manifestations, applications and ways of further development are briefly outlined in Section 5.
2. SPP field parameters based on the microscopic hydrodynamic model.
To start the consideration, we first outline the usual (phenomenological) approach when the materials' constitutive parameters  1 ,  2 =  2 () are postulated prior to the analysis. Then, in each medium, the Maxwell equations hold [24,25] , (here  p is determined via the Drude-model approximation for Ag [23,34]). In this case, the SPP cutoff frequency equals to 1 1 p    = 5.46·10 15 s -1 , and the frequency-dependent behavior of the SPP-field parameters is illustrated by Fig. 2.
The main tasks of this study are the analysis of the solutions (12), (15), (16) and their consequences for the SPP DCs, with special attention to their spatial distributions. For the detailed numerical illustrations, we select the conditions of point A of the dispersion curve in Fig. 2 which corresponds to the radiation frequency  = 4.09·10 15 s -1 and the vacuum wavelength  = 2/k = 0.4609 m. Under these conditions, the metal "effective permeability" 2 = 1 -, SPP wavelength  s = 2/k s and the parameter  become    2 = -2.57,  s = 0.3663 m,  = 6.04·10 7 cm -1 .
(3), (4) and (12), (15), (16). Accordingly, Fig. 3 shows the relative dimensionless amplitudes of the electric and magnetic fields, electrons' velocity and density, normalized with respect to A = E x (0).  (13), (14) and (5)  All curves in Fig. 3a (except the brown one for en/  ) contain large segments of smooth exponential off-interface decay determined by factors exp( 2 x) (x < 0) and exp(- 1 x) (x > 0), which in Fig. 3a look as straight lines. These smooth segments characterize spatial distributions dictated by the phenomenological equations (3), (4); from now on, such contributions and the corresponding terms in the equations will be referred to as the "volume" ones. In the phenomenological description, the "volume" behavior is "in charge" up to the interface (see the blue dashed line for E x as an illustration). But due to the microscopic refinements, just near the interface, the curves' shapes change rapidly, which is caused by the near-surface (NS) terms of Eqs. (12) - (15) proportional to exp(x). The corresponding decay length  is impressively smaller than the "volume" counterparts:    = 0.165 nm, 1/ 2 = 36 nm, 1/ 1 = 96 nm. (19) Note that all the field components are continuous at x = 0 in agreement with the conditions (11) and the 1 st Eq. (17).
The NS contributions are absent in curves for H y and negligible for E z and v z (due to  in the denominators of corresponding terms in (12) and (15)), but the rapid NS variations explicitly enable the continuity of E x and zero value of v x at the metal surface required by the boundary conditions (11). An especially important fact is that the transverse electric field E x vanishes at a certain point x e near the interface (see Fig. 3b; under the conditions (17), (18), x e = -0.21 nm). This fact induces additional peculiar features in the spatial behavior of various field characteristics, as will be shown below (see Figs. 4, 6, 7). In agreement with (16), the curve for the SPP-induced charge density en consists completely of the NS segment (the "volume" value of the non-equilibrium electron concentration is zero). Approaching the interface, the concentration amplitude grows so rapidly that the brown curve goes beyond the scale of Fig. 3a and cannot be seen in Fig. 3b; nevertheless, the condition 0 n n  accepted in the model of (9), (10) is normally not broken: the relative amplitude of the chargedensity oscillations at x = 0 is This equals to 0.739 if A is measured in statV/cm, or 2.22·10 -4 A if A is expressed in V/cm; that is, even in rather high fields, n n (further discussion see in the concluding Section 5). 6 10 A   0 (0) / 1  For the SPP-induced charge density (16), it is important to know its total (or integral) value obtained by the integration over the whole x-range. Such integral values (per unit y-width and unit z-length) make sense for many characteristics, and the unified notation, will be used for such quantities throughout the paper (if necessary, the data for x > 0 can be taken from Eqs. (3) and from the Section S1 of the Supplemental Document). In application to Eq. (16), n(x > 0) = 0 by the definition, and the operation (21) means Now we are in a position to classify the NS terms with respect to their influence on the SPP field characteristics. Eqs. (19) and the curves of Fig. 3a show that, as a rule, the terms of Eqs. (12), (15), (16) with  in denominators ( 1   -terms) make a minor influence and can be safely discarded when compared with the other summands (in the graphs, their contributions are almost undistinguishable). This is not an exclusive consequence of our choice of the medium 2 and of the "working" frequency (point A in Fig. 2) -the Fermi-velocity data for other metals [35] confirm that, in the wide frequency range, the conditions , Nevertheless, the  -1 -terms may produce meaningful contributions upon differentiation with respect to x; that is why their explicit presence in Eqs. (12) - (16) and in some equations below is crucial for the consistent analysis. On the other hand, the role of NS terms can be "suppressed" by the integration, and, generally, the terms proportional to   exp x  make negligible contributions to the integral values (21), except those containing  as a multiplier. As far as (22) holds, such terms can be represented via the delta function: For example, the electron density (16) can be expressed in the form which is actually the limiting form of the electron density dictated by the Drude model [29]. Terms in the form (23), (24) frequently occur in the hydrodynamic description of the medium 2 characteristics when a certain quantity almost vanishes in the volume, but shows an impetuous growth approaching the interface (beside the brown curve in Fig. 3a, see the magenta and green curves in Fig. 5, black and light-green curves in Fig. 7). Due to the delta-shaped behavior, we will call such terms "singular" and use, sometimes, the simplified notation through the delta-function; however, keeping in mind that the boundary values of such terms remain finite, and the deltafunction at any moment can be replaced by (23) with the definite peak value   0    .

Electric polarization of the medium 2
The results (12) -(16) supply a complete characterization of the electron plasma in medium 2 and, particularly, describe its SPP-induced polarization. In the monochromatic field each electron performs an oscillatory motion whose velocity v is associated with the displacement   i  v from its "neutral" position. This displacement produces an effective dipole with a dipole moment per unit volume and Eqs. (12), (15)  "Transverse" and "longitudinal" polarizabilities xx  , zz  depend on frequency . This dependence will be used below, and a useful expression for the derivative can be inferred from (27): Here primes denote derivatives with respect to , and the relations dictated by Eq. (28), are used. In Eq. (29), the second summand vanishes at x = 0 and rapidly decays at x < 0 due to the exponential multiplier, so it merely supplies a minor correction that can safely be omitted in most cases compatible with the condition (22). As for the "longitudinal" polarization   zz x  , as far as its spatial dependence is negligible, the "volume" approximations zz 3. Dynamical characteristics of the SPP field. Now we proceed to calculations of the DCs for the SPP field in the metal layer. In what to come, we distinguish the "field" contributions marked by superscript "F", and the material contributions denoted by superscript "M". The results for medium 1, obtained phenomenologically (see Section S1), will naturally be treated as "field contributions". Also, we will keep the following classification scheme for the spatially-dependent terms:  "volume" terms proportional to exp(2 2 x) and describing the "smooth" exponential decay in the depth of medium 2;  NS terms proportional to exp( x) describing the rapid off-surface decay and giving negligible contributions upon spatial integration;  "singular" terms proportional to  exp( x) describing the strictly-localized near-surface contributions reducible to the delta-function, see (23);   -1 -terms proportional to negative degrees of . These terms are normally negligible and omitted but sometimes kept in intermediate transformations.

Energy density
In the model of Eqs. (9), (10), the energy density in metal is determined by formula [31,33]: the field contribution is expressed by the first summand, the second and third ones describe the material part. By using Eqs. (12) - (16) and neglecting the  -1 -terms we obtain , see Supplemental Document.
First, compare these results to the phenomenological expression (S4). The "field" and material parts are not resolved in (S4) but the (first) "volume" term of (34) expectedly coincides with the phenomenological data. The NS terms contribute to sharp variations near the interface x = 0, well seen in Fig. 4. In contrast to Fig. 3, the NS terms produce not only rapid growth or decrease but also an additional "kink" near the point x e related to the vanishing amplitude E x (see Fig. 3b). The NS terms make no effect on the integral values determined by the "volume" first terms of (32) - (34) but affect the boundary relations. In particular, the energy density "step" predicted by the phenomenological formula (S5) is illustrated by the vertical segment CD in Figs. 4a, b; the "microscopically corrected" step following from (34), characterized by the segment CE, is higher. Interestingly, the "field" part of the energy density (blue curve) is distributed continuously (due to the continuity of the electric and magnetic fields, see (11) and (17)), so the complete energy discontinuity (35) follows from the 2 nd and 3 rd material terms of Eq. (31), which exists only in the medium 2. Now consider the energy flow density. In the microscopic approach and in application to the TM SPP field, it is given by equation [33] It contains the material part additional to the field contribution (S7) and expressed by the second summand of (36). According to (13) - (16), this material part can be presented as  which fully consists of γ -1 -terms (see the classification of the 1 st paragraph of Section 3), and, due to our convention, can be discarded, whereas the first summand of (36) via the immediate application of second Eq. (4) and Eq. (12) leads to the expression The first term of the right-hand side coincides with the phenomenological result (S8) but due to the NS summand of (37), the energy flow step at the interface differs from (S9): In the example of Fig. 4, the "phenomenological" step AB is reduced to because the electric and magnetic fields are continuous ( 1 = 1, see (17) and (38)). Besides, Eq. (37) shows that the energy flow inversion occurs inside the medium 2 rather than exactly at the interface. The "inversion point" is the same point x e where E x changes its sign (cf. Figs. 3b and 4b) -quite expectedly in view of the definition (S7), (36). And the general behavior of the microscopic energy flow in the NS region looks rather similar to that of the transverse electric field component. The integral energy flow is still determined by (S10).  (34) and (37)). D and A: surface values of the phenomenological distributions (S4), (S8) for x < 0; C and B: surface values of distributions (S4) and (S8) for x > 0. The interface plane is marked by the vertical dashed line, the panel (b) shows the NS region of (a) in a magnified horizontal scale.

Momentum density.
Our next step is the microscopic calculation of the momentum in medium 2. In contrast to the phenomenological analysis (S16) -(S24) [23], which was based on the standard formal procedures, the microscopic approach resembles the process of "building" the "complete" momentum from certain "blocks" of different physical origins.
The first momentum constituent is evidently the "pure-field" contribution expressed by the Poynting vector In view of (S7), this quantity is closely related to the energy flow density (37), In principle, its behavior can be understood from the above discussion of S ( Fig. 4) but in the "momentum context" it is again illustrated by the blue curves in Fig. 5. Remarkably, the "volume" part of (40) perfectly coincides with the volume part of the phenomenological momentum (S23).
Other "blocks" of the electromagnetic momentum in medium 2 include the material contributions associated with the motion of electrons. To find these contributions, we must consider the force exerted by the field on the medium particles, which is performed in the Supplementary section S2 and dictates that the material contribution consists of two parts. The first one is determined by Eq. (S29) and can be presented as where, in view of Eqs. (13), (16) and (4), and F p is determined by (39), (40) so that finally Its spatial behavior is visualized by the yellow curves in Fig. 5. Note that the NS term causes it to vanish at the boundary. "Block" (41) - (43) represents the material part of the Poynting momentum proportional to  = 1 - 2 [16].
The second material "block" of the SPP momentum, which can be called "electric" due to direct connection to the electric field components, follows from (S30) in the form For the SPP field of Eqs. (4), (12) and employing only the 1 st term in the right-hand side of Eq. (29), this can be expressed as This approximation is shown by the orange curves in Fig. 5. For comparison, the results of a direct application of Eq. (44) with numerical differentiation for xx  and zz  are presented by the cyan curves: the difference may reach ~10% but it is localized in the very narrow NS region and practically does not affect the general character of the spatial dependence (45). Yet another material contribution appears due to the rotational motion of electrons. Its calculation starts from the "material rotational spin" [20][21][22]31]: which, via Eq. (15), entails Hence, the "material rotational spin momentum" follows from the spin momentum definition (S20) in the form M where the relation (23) has been used; note that the 1   -terms are held in (47) because these make meaningful contributions in (48). The expected boundary value of the "volume" part of M SR p (given by the first summand in brackets of (48)) is p/|A| 2 , 10 -12 s/cm (see point A in Fig. 5a). The "regular" NS contributions (2 nd and 3 rd terms in brackets of (48)) are practically unperceivable on the background of the "singular" delta term; the latter appears due to the fast decay of the spin (47) from its volume behavior (1 st summand in parentheses (47)) to zero, caused by the second (NS) summand (see Fig. 8 which is well illustrated by the green curve in Fig. 5a. As a result, the overall electromagnetic momentum p in the medium 2 (not shown in Fig. 5) appears as a sum of the pure-field contribution (40) (blue curve in Fig. 5) and the "singular" momentum (49) with a near-vertical ramp of the green curve in the interface vicinity: Noteworthy, both the "volume" and "singular" terms perfectly agree with the phenomenological expression (S23), and only the NS terms allow for refinements inspired by the boundary conditions for E x and v x (11) in the microscopic theory. This agreement confirms that all "blocks" (40), (43), (45), and (48), which form the SPP momentum in medium 2, are identified correctly, and the microscopic theory complies with the phenomenological approach presented in Supplemental Document. Additionally, the microscopic analysis discloses the nature of the "singular" spin momentum surf p of (S21) and shows that its "true" boundary value is finite: [which equals to 5.05·10 -10 g/(cm 2 ·s) with conditions (17), (18)]. In this expression, non-"singular" terms of (48) -(50) are neglected because, due to (19), the right-hand side of Eq. (51) exceeds the other momentum constituents by two orders of magnitude. Because of the high absolute values, only the initial stages of the exponential "ramp" of the green and magenta curves, approaching the interface, are shown in Fig. 5a.

Spin-orbital momentum decomposition
In the phenomenological approach of Section S1 (see Supplemental Document), the spin-orbital (canonical) momentum decomposition is basic for the momentum calculation. In the course of the microscopic approach it is more natural to rely on another meaningful subdivision of the momentum constituents -namely, separation of the material and "field" contributions, which was made in the previous subsection. Now we consider the spin-orbital decomposition in the microscopic framework. To this end, we separately inspect each constituent of the field momentum calculated in subsection 3.2, We begin with the first terms whose attribution is quite evident. One can easily see that expression (44), upon discarding the NS terms of (29) and taking (30) into account, can be represented as that is, M e p implies the dispersion correction of the orbital momentum (S16) and is therefore its part. On the other hand, M SR p , originating from the spin (46), belongs to the spin momentum "by the definition".
Other terms need a more careful investigation. Regarding the pure-field momentum p F , its decomposition immediately follows from the formal decomposition of the Poynting vector (39), which in view of the Maxwell equations (9) reads [20,22] This form differs from the usual spin-orbital decomposition [27,28,30] valid only for fields without free charges and currents. The presence of charges in Eqs. (9) destroys the "electric -magnetic democracy" [27] and produces additional charge-dependent terms (last summands of (54) and (55)) whose attribution is not dictated by general principles and can be ambiguous. The form (54), (55) is heuristic and apparently contradicts to some natural expectations: e.g., the expression (54) leads to a non-zero integral value and thus violates the spin-momentum property (S22). Such controversies emerge because in presence of free charges the Poynting momentum (39) is only a part of the entire electromagnetic momentum. As we will see, they are removed after involvement of other momentum contributions, and the decomposition (53) which show that the "singular" -proportional terms of both summands mutually cancel. Since the delta-shaped ( -proportional) "singularity" of the "field" spin momentum at x = 0 is hardly interpretable and looks non-physical, such a cancellation testifies that the term ~  * Im n E has been reasonably attributed to the spin part of the decomposition. As a result, an explicit expression for the spin part of the "field" momentum (54) follows in the form In ( Note that the "orbital-like" form of this term (cf. Eq. (S16)) confirms its attribution to the orbital part (55) of the decomposition (53). As a result, the "field" orbital momentum (55) acquires the form By using Eqs. (5), one can easily verify that the sum of (56) and (57) coincides with (40), and that the "volume" terms coincide with the phenomenological results (S18), (S21). At last, the decomposition of M m p can be performed quite similarly due to relation (42) connecting it to the "field" momentum F p . The only problem that remains is to separate properly the spin and orbital parts of the additional NS summand in (42), which is solved by the natural requirement for both parts to vanish at the interface. Accordingly, M Interrelations between the elements of the spin-orbital decomposition and their spatial distributions are illustrated in Fig. 6 where they are confronted with the corresponding momentum constituents in medium 1 determined by (S18), (S21). One can see that the "volume" orbital and spin constituents are directed oppositely, so the complete momentum appears as a result of their "competition". Notably, in medium 1 the orbital momentum F O p is "stronger"; its direction coincides with the complete momentum direction but the absolute value is higher,   6 shows that while the phenomenological results ("volume" terms) stipulate a behavior discontinuous at the interface for all the considered quantities, the NS terms restore the continuity of the momentum constituents. Indeed, the blue, green and red curves of the "field" contributions F O p , F S p and F p are continuous with their counterparts in medium 1 (herewith, the continuity of F p is "absolute" and follows from the definition (39) and the boundary conditions for E and H (11), whereas F O p and F S p are continuous only in the considered situation  1 = 1). The material spin and orbital constituents vanish at x = 0, which can also be treated as a sort of continuity because the material contributions in medium 1 are absent. The data of Fig. 6 illustrate mainly the technical details of the spin-orbital decomposition. In addition, Fig. 7  which include not only the elements of the formal spin-orbital decomposition (53) -(55) but also the terms M e p and M SR p (see (52), (45), (48)) that were classified at the beginning of this Section based on their immediate physical meaning. Then, by combining the "field" and material contributions, we build the complete spin and orbital momenta of the field, It is expressions (62) and (63) that should be immediately juxtaposed with the phenomenological results of Eqs. (S18) and (S21) but even Eqs. (60), (61) 'per se' permit to draw some physically consistent conclusions. In particular, they tell that the "volume" parts of M O p and M S p have the same absolute values but are directed oppositely (the "volume" parts of the light-green and cyan curves in Fig. 7a are symmetric with respect to the horizontal axis). That is why in the complete material momentum p M they exactly compensate each other, due to which Eq. (49) only contains the NS and "singular" terms. Remarkably, despite that the "volume" contributions of M O p and M S p mutually cancel, the integral material contributions to the orbital and spin momenta over the whole medium 2 are equal: (… 2 means the integration like (21) but limited to the x < 0 half-space). This happens due to the "singular" term in M S p (61) that "acts" oppositely to the "volume" one and produces a twice higher integral contribution.
Further, cyan and magenta curves in Fig. 7a demonstrate, in full agreement with Eqs. (59) and (60), that the "volume" parts of M O p and M mS p coincide, but in the NS region these quantities behave quite differently (see Fig. 7b). Generally, in the NS region all the curves show characteristic "kinks" originating from the transverse-electric-field sign inversion and whose nature was already discussed in connection to Figs. 4 -6. Exclusions are the black and light-green curves presenting the material M S p and the complete p S spin momenta. In these dependencies, the dominating role belongs to the "singular" -shaped term (originating from the material "rotational" spin momentum M SR p , see the magenta curve in Fig. 5a) whose sharp growth masks all fine details. To make them visible, the "residual" spin momentum  which is the complete spin momentum with exclusion of the "material rotational" constituent (48), has been calculated, and the result is illustrated by the yellow curve in Fig. 7, where the "kink" is distinctly seen.  (17), (18): Spin-orbital decomposition of the "field" (40) and material (43) momentum. (Blue) complete "field" momentum (40) and (S23) (for x > 0); (green) orbital "field" momentum (57) and (S18) (for x > 0); (red) spin "field" momentum (56) and (S21) (for x > 0); (cyan) orbital material contribution (58); (magenta) spin material contribution (59); (yellow) material Poynting contribution (43). Blue and yellow curves are the same as in Fig. 5. The interface plane is marked by the vertical dashed line, the panel (b) shows the NS region of (a) in a magnified horizontal scale.
Note that only the illustrative, artificial quantity p at x = 0 (see Fig. 6). All other momentum constituents discussed in Fig. 7, including the complete orbital and spin momenta p O (62) and p S (63) experience "jumps" at x = 0, which can be attributed to the fact that these are formed as combinations of the "field" and material parts, whereas in medium 1 only the "field" contributions exist. Noteworthy, the NS terms destroy the perfect proportionality between the orbital momentum p O (62) and energy density w (34) distributions that are characteristic for the phenomenological approach [23] (cf. Eqs. (S17) and (S4), (S18)): their difference in the NS region is shown by the grey and brown curves in Fig. 7b. However, the boundary values p O (0) and w(0) in Fig. 7b visually coincide (this is an artifact of the supposition  1 = 1 (17)).
Finally, one can see that the "volume" parts of p S and p O perfectly agree with the phenomenological expressions (S18) and (S21) as well as the "singular" term of (63) agrees with surf p in (S21). Accordingly, all integral (in the sense of (21)) properties of the phenomenological expressions hold for the microscopic model. In particular, the microscopically defined orbital momentum obeys integral relation (S19), and the integral spin momentum is equal to zero (as it is dictated by Eq. (S22) and required by the general theory [27,28]) despite that it was derived with the use of intermediate spin-momentum constituents based on Eq. (54) which apparently violates this requirement. As to the "surface" -term of (S21), the microscopic analysis reveals its "true" finite boundary value (see (51)) which, due to (19), is two orders higher than other momentum contributions described by Figs. 5 -7 and is thus not shown in the figures. Again discarding the relatively small non-"singular" terms as it was made in (51), we find that the result (51) is also correct for . Fig. 7. SPP momentum constituents calculated for the conditions (17), (18): Complete spin and orbital contributions and their material ingredients. (Blue) complete "field" momentum (S23) (for x > 0); (green) orbital "field" momentum (S18) (for x > 0); (red) spin "field" momentum (S21) (for x > 0); (brown) complete orbital momentum (62); (black) complete spin momentum (63); (cyan) orbital material momenum (60); (light-green) spin material momentum (61); (magenta) spin material contribution (59) (the same as in Fig. 6); (yellow) "residual" spin momentum (65); (grey) energy density (34) multiplied by (k s /) according to (S17). The interface plane is marked by the vertical dashed line, the panel (b) shows the NS region of (a) in a magnified horizontal scale.

Spin density of the SPP
A natural way for the microscopic spin calculation is formally putting together different contributions associated with the rotational motion of the field vectors and field-driven material  [31,20,22]. Formally, this way requires a detailed inventory of the field-induced rotations in the medium, similarly to what was performed upon considering the momentum contributions in subsection 3.2. One of the spin contributions has already been mentioned: it is the material spin associated with the rotational motion of electrons (46), (47), which was employed as an auxiliary item for the momentum analysis. Here we reproduce this contribution in the final form with omitted 1   -terms: This "rotational" material spin is illustrated by the cyan curves in Fig. 8. It satisfies the natural condition : the very rapid variation from the expected "volume" value A to zero is caused by the NS term. This term is also responsible for the "singular" surface component of the momentum (48). Another spin contribution is the "field" spin which formally follows from Eq. (S12) applied to the microscopic TM field in medium 2 with  =  = 1 and obeying Eqs. (9). However, the "formal" spin definition (S12) was obtained for the fields without free charges, and its non-critical use may miss some medium-associated corrections that were already applied for the momentum calculations above. To avoid this, the corresponding spin contributions will be calculated based on the relation (S20) between the spin and spin momentum due to which, for the SPP geometry, any z-oriented spin momentum "block" p S is related to the associated y-directed spin contribution s: cf. the blue and green curves in Fig. 8. It is instructive to compare this result with the "formal spin" derived from the formal spin definition (S12) (light-green curves in Fig. 8). First to note, reduces to the "field" spin (however, modified by the NS term): the specific material contribution  (Fig. 8) while the spin rapidly changes near the interface, reversing the sign together with E x in point x e , cf. Fig. 3b. The NS term in (69), responsible for this difference, appears because the term ~F ormal s   * n F s Im of (54), taken into account in (56) and, consequently, in (68), is omitted in (69). As a result, is discontinuous at x = 0 while is continuous with the spin in dielectric (Fig. 8) -in full compliance with the boundary conditions (11) for the electric field. However, another case of continuity is more important: the quantity is just the "naive" Minkowski spin (described by Eq. (S13) at 2 2     ). It satisfies the condition (see the red curves in Fig. 8), which agrees with the fact that the Minkowski spin is continuous [23] (cf. also Eq. (S13) and the note beneath (S14)). This conclusion supplies additional arguments to the correctness of Eqs. (68), and justifies the procedure chosen for the spin calculation via (56), (59) and (67) -had the spin been derived formally as (69), the Minkowski spin continuity would have been destroyed.
The "complete" spin in medium 2 (see the magenta curves in Fig. 8) is described by the equation Its "volume" part coincides with the phenomenological result (S13); however, the NS term in (71) provides that, in contrast to the phenomenological spin of Eq. (S13), the overall microscopic spin is continuous: This formula replaces the phenomenological result (S14) depicted by the "step" BC: the microscopic analysis enables both the "naïve" Minkowski spin and the complete spin to be continuous at x = 0. Another interesting observation concerns the "spin-momentum law" [36] which interprets the "extraordinary" transverse spin of the SPP through the vorticity of the longitudinal momentum density [6,36]: (see Eq. (S25)). This feature of the SPP field reveals a transparent physical analogy with the layered flow of a fluid near a plane wall, which is accompanied by the flow vorticity orthogonal to the main direction of propagation [37]. However, juxtaposition of Eqs. (71) and (50) indicates that the "spin-momentum law" (S25) is only valid for the volume contributions, while close to the interface the proportionality between the spin s and the curl of momentum p is destroyed by the NS terms.
   s p 4. Frequency dependence of the SPP dynamical characteristics The results of Sections 2 and 3 provide an exhaustive description of the spatial distributions for various DCs of the SPP field. However, the illustrations were restricted to the single radiation wavelength corresponding to point A of the dispersion curve (Fig. 2). In this Section, we briefly outline the main modifications expected when the radiation frequency  or, which is the same, the vacuum wavelength , changes.
We start with the remark that, according to Fig. 2, for any frequency the parameter  (13) remains two orders of magnitude higher than other spatial parameters  1 and  2 , and all approximations associated with the condition (22) are thus valid within the wide range of possible frequencies. Consequently, all the analytical formulas of Sections 2 and 3 are appropriate, and the frequency dependence of the DCs characteristics is comprehensively determined by the frequency dependencies of the main SPP parameters k s ,  1 ,  2 presented in Fig. 2. Fig. 9 provides a comparison between the SPP characteristics for the two situations: "highfrequency" regime specified by point B in the dispersion curve of Fig. 2 ( = 5.19·10 15 s -1 ,  = 363 nm) and the "low-frequency" one (point C in Fig. 2,  = 2.685·10 15 s -1 ,  = 702 nm). The first notable discrepancy between the two situations is remarkably higher absolute values of all the presented quantities in the high-frequency regime. However, this observation is not new: As it is shown in previous sections, the "volume" contributions as well as the integral values (in the sense of (21)) are correctly described by the phenomenological expressions of Section S1 (see Supplemental Document) whose frequency-dependent behavior is well known [1,2,23]. Therefore, here we focus on the NS contributions inducing some peculiar features of the DCs' spatial distributions.
For many DCs (e.g., those discussed in Figs. 5,8), the frequency variation entails mainly quantitative transformations of the NS behavior coherent with changes in the absolute values of the "volume" terms. But more interesting NS features are associated with the transverse-electric-field sign inversion (which is rather frequency-sensitive, see the "high-frequency" x eB and "lowfrequency" x eC zero-point positions in Fig. 9a, b and Fig. 10a). These are the characteristic "kinks", additional "peaks" and "troughs" in the distributions of energy and of the elements (58) -(61) of the material-momentum spin-orbital decomposition (see Figs. 4, 6, 7). The corresponding momentum constituents are presented in Fig. 9, and one can see that the shapes of these peculiar details substantially depend on the radiation frequency. Fig. 9a shows that the "potential wells" in the NS energy distributions, distinctly seen at the "high" (and "moderate", see Fig. 4) frequencies, disappear when the frequency decreases. Interestingly, this happens because their "left walls" w 2 , fall down whereas the "right walls" w 1 , practically do not change. In more detail, the frequency dependence of the "left walls" is illustrated by Fig. 10a for the "field" and "complete" energy density (the corresponding data for the can be derived from Eq. (31)).   x, nm Similarly, Fig. 9b juxtaposes the "low-frequency" and "high-frequency" "peaks" and "troughs" of p O , M O p and M mS p distributions (cf. Fig. 7). In contrast to Fig. 9a, where the "well" centers (energy minima) move "inwards" the medium 2 together with x e (cf. the black curve in Fig. 10a), the extremum positions of the "kinks" are less frequency-sensitive and always lie closer to the interface. However, the "kinks" morphology characterized by the parameters p O , δ M O p and δ M mS p explained in Fig. 9b shows more interesting frequency-dependent evolutions presented in Fig. 10b.
The data of Figs. 9 and 10 supply a general characterization of the frequency-induced changes in the NS spatial distributions of the selected DCs. In particular, they show possible ways for the deliberate control of these spatial distributions, especially, for purposeful formation of the "kinks", "peaks" and "troughs" of the corresponding curves inside the nanometer-scale vicinity of the interface.

Concluding remarks
The main outcome of this paper is the detailed analysis of the near-interface spatial behavior of the DCs in a lossless SPP wave. Due to the use of the hydrodynamic model, extremely subwavelength structures are revealed in the distributions of the field energy and momentum. Importantly, the paper presents the first to our knowledge collection of analytical expressions describing the separate constituents of the electromagnetic momentum in plasmonic fields, with the explicit identification of the "field" and material contributions of different origins. For a consistent comparison with the known phenomenological results, which are expressed in terms of the spin-orbital (canonical) momentum decomposition, the procedure of such decomposition in presence of free charges and currents is developed. The procedure validity is substantiated by using the SPP fields as a model example.
Within the limits of common applicability, the new details of the DC spatial distributions comply with the known SPP properties described by the phenomenological approach and/or by the simplified Drude model of the metallic medium. However, the new results demonstrate the existence of controllable nano-scale optical structures in the NS layer of the medium 2 that can be deliberately created and whose parameters can be efficiently regulated, e.g., by choosing the SPP frequency.
Although the present consideration deals with the homogeneous metallic medium 2, it can be generalized to more complex plasmonic structures containing the metamaterials, composite media, fluids, nanoclusters, etc. [39][40][41][42]. One can expect that specific local interactions between the field and material nano-objects (e.g., impurity atoms, molecules, lattice defects and nanoclusters introduced into the metal matrix [43]) can be inspired and mediated by the nanoscale field inhomogeneities described in the paper, and the degree of spatial commensurability between the field "kinks" and the nano-objects supplies additional channel for controlling this interaction. The subwavelength field structures can be useful in various problems of nano-optics: for example, the "peaks" and "troughs" in the NS energy distribution (Figs. 4, 9a) can serve for precise localization of the objects, their spatially-selective excitation or transformation, etc. Usually, such nano-objects are characterized by a noticeable mobility in the crystal structure, or the medium 2 can be liquid [44]; then, the controllable "kinks" in the momentum distributions (Figs. 6, 7, 9b) may be employed for their guiding and transportation due to the field's mechanical action [38].
Especially interesting manifestations of the specific DC-distribution details predicted in this paper are expected in combination with the external or SPP-induced static fields, for example, the SPP-induced magnetization [20][21][22][23], which is proportional to the material spin constituent M R s (66). In particular, the nanoscale variations of the energy and momentum distributions may affect the electron spin transport processes [45][46][47], induce additional influences on the photo-sensitive centers in presence of high-gradient local optical fields [48], etc. Despite the presence of a non-zero spin density, which is a hallmark of SPP fields [7,8,11,12], the usual TM SPPs produce no specific chiral action because their transverse spin is helicity-independent. However, the use of specific chiral metasurfaces [49,50] enables generation of chiral SPPs with controllable spin orientation, which can selectively excite and control chiral nano-objects, e.g. the diamond spin-qubits [51].
In conclusion, let us outline the main restrictions of the presented approach and possible ways of its further development. The first limitation originates from the lossless approximation, which is widely used in the SPP theory [2][3][4][5][6][7][8]12,13,20,23] and enables the physically transparent analytical description. In principle, the effects of damping (energy dissipation) can be taken into account by an additional summand in the right-hand side of Eq. (10) where  is the collision frequency [30,32]. After this modification, the formal solution of Eqs. (1), (9), (10) and (11) which, in turn, lead to the modified expressions (7) and (13) As a result, the expression for  2 (7) is transformed to the usual complex permittivity expression for the Drude model with losses [2], and the constant  becomes complex. Unfortunately, such direct modifications cannot be applied to the DC expressions because these are quadratic in field and electron velocity components and include the real and imaginary parts separately; however, some general conclusions can be inferred from the simple reasoning. In order the propagating SPP to exist, the damping correction must be small; for example, in the considered case of silver,  = 5.9·10 13 s -1 [52] so for the conditions of points B, A and C of Fig. 2, / = 0.049, 0.014 and 0.022, respectively. Such small corrections mean that in the most interesting cases, the loss-induced modifications of the results are mainly quantitative. Nevertheless, some corrections are principal.
For example, due to complex  and  (72), the electric field component E x (12) differs from zero everywhere in the medium 2: while Re(E x ) vanishes at x = x e , the small "imaginary" field E x , which is nearly in-phase with the longitudinal component E z , still exists. Accordingly, the singularity in xx  (27) is removed. Practically this means that the "kinks" in the energy and momentum distributions become smoothened, as well as the "singular" NS terms in Eqs. (16), (48) - (50). Anyway, the detailed analysis of the DC behavior in presence of energy losses deserves a special investigation in a separate work.
The second important simplification is the assumption of perfectly flat interface, which is frequently broken in practice [53][54][55]. Usually the influences of roughness are evaluated by statistical methods but, for the effects considered here, this way of reasoning cannot be applied because even very small boundary imperfections of the order of |x e | (see Figs. 3 -8) are crucial. Most probably, the irregularities of such size scale will induce local SPP-field confinement and enhancements [29] but their correct analysis can hardly be performed analytically and requires exquisite numerical approaches (e.g. the finite-element method with non-uniform meshing [29]). One may expect that the out-of-plane boundary deviations of the size ~|x e | will not destroy the predicted pictures as far as their size in the (y-z)-plane ( Fig. 1) exceeds |x e | by the order of magnitude. Alternatively, the interface irregularities can be used for purposeful formation of "hot spots" in the SPP fields [53].
The concept of the "interface imperfections" makes a natural "bridge" to the complex SPPsupporting structures where the media 1 and 2 are separated by the specially created metasurfaces, contributing, for example, to the formation of chiral SPPs [49,50]. Usually, such meta-boundaries form the physically selected nanometer-wide intermediate layer between the metal and dielectric media, which in calculations is replaced by the conventional plane interface with specific electromagnetic properties [50]. The applicability of the existing hydrodynamic model in the vicinity of this conventional interface is questionable but its further generalization to such situations looks promising and relevant.
And the last approximation, which deserves a few remarks, is the linearity of the basic equations (1) and (9) - (11). In contrast to the usual SPP theory where the linear character is accepted as needless to say, in this paper the situations of extraordinary high local values of some parameters occurs systematically, e.g., in connection with the "singular" NS terms in Eqs. (16), (48) - (50). This requires a special attention: for example, if in Eq. (20) the non-equilibrium concentration of electrons n(0) exceeds or even becomes comparable with the "background" concentration n 0 , this is physically inappropriate. Estimations made in Sec. 2.2 show that the requirements to the admissible field intensities are rather loose. If we consider an SPP wave with the total power Q propagating within the width y, the electric field amplitude can be determined via Eq. (S10) as    . This value, probably, warrants the model applicability in many practical situations; however, in higher fields which can easily occur due to microscopic dimensions of the system [29], the improvements accounting for the quantum character of the electrons' behavior, including the quantum tunneling effects [29] as well as the non-linear response of the media [4,56], would be necessary. Microscopic analysis of the energy, momentum and spin distributions in a surface plasmon-polariton wave Supplemental Document S1. Phenomenological characterization of the SPP field S1.1. Energy distribution.
Taking into account dispersion of the medium, the energy density is determined by the Brillouin's formula [24,25]