Surface scattering amplitude for a spatially dispersive model dielectric

Recently, it was shown [H.-Y. Deng and E. A. Muljarov, Phys. Rev. B 106 , 195301 (2022)] that the so-called additional boundary conditions (ABCs) had signiﬁcant issues when describing the surface properties of spatially dispersive media and an alternative ABC-free theory was developed introducing the surface scattering amplitude (SSA), which determines how polarization waves are scattered by a surface. Here we analytically calculate the SSA for a spatially dispersive model dielectric using wave-scattering theory and discuss some emerging generic properties independent of the model. The model consists of a lattice of interacting harmonic oscillators which had historically been used to mimic local excitons. It was shown that for short-range interactions the SSA is a constant independent of where the polarization waves were originated, which is not the case for long-range interaction. The calculated optical properties of a slab verify the general ABC-free theory.


I. INTRODUCTION
It is a long known fundamental issue in macroscopic electromagnetism that Maxwell's boundary conditions (MBCs) are incomplete for determining the optical responses of spatially dispersive media [1][2][3].In such a medium, the dispersion relation admits more than one optical mode of the same frequency and polarization and the MBCs are insufficient to uniquely determine the optical reflection and transmission amplitudes.The prevailing remedy, first proposed by Pekar in 1957 [4,5], has been to supplement MBCs with additional boundary conditions dubbed ABCs [6][7][8][9].However, the justification of ABCs remains an interesting topic up for debate [10][11][12][13][14][15][16][17][18][19][20][21][22][23][24][25][26][27].Recently [28], it was shown that ABCs furnished incorrect descriptions of the physical attributes of the surfaces of spatially dispersive media.For example, ABCs cannot be surface characteristics and fail to reproduce size quantization effects.The correct description was developed introducing surface scattering amplitude (SSA), which serves as the fingerprint of a physical surface.The SSA for a surface describes how incident polarization waves excited by a perturbation are scattered by the surface.It can be both experimentally measured and theoretically calculated [28].
Here we calculate the SSA for a spatially dispersive model dielectric and study its properties.The employed model has widely been used to describe, for example, the optical response of Frenkel excitons [15].It consists of a cubic lattice (with lattice constant a set to be unity throughout this work) of harmonic oscillators, one at each site labeled by n = (n x , n y , n z ) with integers n α , where α = x, y, z.The displacement of the oscillator at site n is denoted by r(n, t ) = [r x (n, t ), r y (n, t ), r z (n, t )].An oscillator is assumed to carry an electric dipole d(n, t ) = er(n, t ), where e is an electric charge.Polarization waves (i.e., dipole density waves) can be excited in the presence of light.The classical equation of motion for the oscillators reads [15] L αβ (n, m)r β (m, t ) + e M E α (n, t ). ( Here α, β = x, y, z, ω 0 is the natural frequency of the oscillators, M their effective mass, and γ their friction coefficient, E α (n, t ) denotes the αth component of the electric field E(n, t ) present at site n, the quantity L αβ (n, m) gives the coupling between oscillators, and (n) brings about irregularities near the boundaries of the lattice and describes the boundaries themselves.We assume homogeneous and isotropic coupling and hence L αβ (n, m) = δ α,β L(n − m), where δ α,β denotes the Kronecker symbol.By the value of , a bounded system naturally divides itself into one or more surface regions (SRs) and a bulk region (BR), where the SRs comprise all sites with = 0, while the BR includes all the rest [28].One may create a bounded system by rendering in an unbounded (infinite) system a set of sites with = ∞, which effectively isolates a portion (e.g., a slab) from the rest.These virtual sites may also be included in the SRs.Note that on these virtual sites r = 0 as the oscillators possess infinite stiffness.
In the present work we are concerned with systems possessing planar boundaries (i.e., flat surfaces), i.e., semi-infinite systems and slabs.By virtue of the (discrete) translational symmetry along the surfaces which are assumed parallel to the x-y plane, it suffices to consider plane waves, where n = (n x , n y ), k is the planar wave vector, and ω is the frequency.Hereafter we shall write n z as n to avoid notation cluttering, where n labels a layer of sites (see Fig. 1).Assuming that only depends on n, we rewrite Eq. ( 1) as where ω2 = ω 2 + iγ ω and and hence we ignore it subsequently.
It suffices to focus on one component so we take and d we can rewrite Eq. ( 2) as In the rest of this paper, we assume that L k (n − m) vanishes for |n − m| > l c , where l c represents the coupling length.For l c = 1 there is only nearest-neighbor coupling.
In the next section, we solve the equation of motion and obtain the SSA.Then the optical properties of the system are studied in Sec.III and the conclusion is presented in Sec.IV.In Appendix A we provide a proof that virtual sites cannot be excited and can be used to create boundaries.In Appendixes B and C, analytical results are reported for hard-wall surfaces with nearest coupling (l c = 1).

II. SURFACE SCATTERING AMPLITUDE
A polarization wave of unit amplitude incidents upon a SR and gets reflected with an amplitude which we call the SSA [28].In this section, we calculate the SSA for the model medium described above by solving the equation of motion (2) or equivalently (4).We shall do this by first solving the equation for an infinite system for which H = H 0 and then for a semi-infinite system for which only one SR is present and finally proceed to solve that for a slab possessing two SRs.A schematic of the slab system is shown in Fig. 1.We create such systems by rendering a stack or stacks of virtual layers which decouple the system of interest from the rest; see below for details.
Before solving the equation, it is worth noting that, for an infinite system, Bloch's theorem dictates that the solution in the absence of the electric field (i.e., for E = 0) represents waves of exponential form, i.e., r(n) ∼ e iqn , where q satisfies the following equation: This equation determines the dispersion of polarization waves existing in the system.The number of total possible values of q for a given ω equals 2l c , since Eq. ( 5) is a 2l c th order polynomial equation of Z = e −iq and thus has 2l c roots.The resulting q ν are in general complex.With parity symmetry, i.e., L k (l ) = L k (−l ) assumed throughout the paper, one-half of these roots have positive imaginary parts and correspond to right-going waves (i.e., waves propagating or decaying to the right) and the other half (with negative imaginary parts) correspond to left-going waves.We denote by q ν those for right-going waves, where ν = 1, . . ., l c , thus leaving −q ν for left-going waves.

A. Infinite system
For an infinite system without SRs, we allow n to go from −∞ to ∞ and have H = H 0 .The equation of motion becomes We neglect self-sustained waves (i.e., the solutions to the homogeneous part of the equation) and formally solve this equation as where S 0 (n, m) gives the wave amplitude r(n) if (e/M )E (n) = δ n,m and is the Green's function of H 0 , which satisfies Requiring that S 0 describes out-going waves and making use of parity symmetry, one can easily show that [29] S 0 (n, m) = where S ν are coefficients independent of n and m, which are determined by applying Eq. ( 8) to n = m + 1, . . ., m + l c .For nearest-neighbor coupling, l c = 1 and there is only one mode.The corresponding wave number we simply denote as q and the coefficients as s.We find S 0 (n, m) = s e iq|n−m| , s = 1 2iL 1 sin q , L 1 = L k (1), (10) with ω2 = ω 2 0 + 2L 1 cos q.

B. Semi-infinite system
Formally, a stack of virtual layers can divide an infinite system into two semi-infinite subsystems, which are isolated from each other as long as the stack is thick enough (i.e., consisting of at least l c virtual layers; see Appendix A).For the sake of definiteness, we choose a stack consisting of l c virtual layers.Let us then focus on one such subsystem that extends from n = 1 − l c to n = ∞, including the virtual layers (i.e., 1 − l c n < 1).The SR then consists of the layers n = 1 − l c , . . ., j c , where j c is the physical width of the SR, and the remaining layers comprise the BR.The equation of motion is given by Eq. ( 4), with H being that for the semi-infinite system.Again, we look at the special case where (e/M )E (n) = δ n,m and find the corresponding Green's function, which is denoted by S(n, m).It satisfies The general solution to Eq. ( 4) can then be written as where the self-sustained part has been left out so that only excited waves are considered.
To proceed, we resort to a general relation (Lippmann-Schwinger relation [29]) between S and S 0 , that is, S = S 0 + S 0 T S 0 , where the scattering matrix describes the scattering effects due to H .It follows that where T (l, l ) are the elements of the scattering matrix T .As H (l, l ) vanishes for l and l outside the SRs so must T (l, l ).As such, the summations over l and l in Eq. ( 14) can actually be restricted to layers in the SRs.One may note that S is a symmetric matrix if H and (hence T ) is symmetric.In Eq. ( 14), the first term is the same as for an infinite system without any scattering effects, while the second term contains all scattered waves due to the boundaries (i.e., the SRs).The second term describes the process of a wave excited at layer m by an electric field to propagate to layer l in the SRs, where scattering activates layer l to generate secondary waves which propagate to layer n.In Appendix A, we show that S(n, m) ≡ 0 whenever n or m or both correspond to virtual layers, which is expected since virtual layers cannot be excited and, as a result, S(n, m) vanishes also whenever n and m are on opposite sides of the stack of virtual layers.
For the semi-infinite system, S(n, m) takes on a particularly simple form for n in the BR.In such case, n lies to the right of the SR and S 0 (n, l ) = ν S ν e iq ν (n−l ) , where l j c < n is assumed.Hence the second term in Eq. ( 14) represents purely right-going waves ∼e iq ν n , which are nothing but the reflected waves due to the SR.We then obtain Here S ν is the same as for an infinite system defined in Eq. ( 9) and R ν (m) is by definition [28] the surface scattering amplitude (SSA) for the νth mode.The factor e iq ν m in the second term is picked up by a wave (in the νth mode) traveling from the layer m, where it is excited, to the SR.
The as-defined SSAs, R ν with ν = 1, . . ., l c , are the fingerprints of the SR.Substituting Eq. ( 9) into ( 14) leads to This expression admits a simple interpretation [30]: an incidence wave excited at layer m travels to the SR and drives a layer l therein [i.e., the factor S 0 (l , m)], thereby activating layer l [i.e., the factor T (l, l )] to generate a reflected wave with amplitude R ν .The factor e −iq ν m cancels the factor e iq ν m in Eq. ( 15), while the factor e −iq ν l offsets the phase due to the location of the layer l.The SSAs in general depend on m, thus providing an intrinsic mechanism for the SSA spatial dependence introduced in Ref. [28].However, for hard-wall surfaces ( j c = 0) with nearest-neighbor coupling (l c = 1), which is studied in detail in Appendixes B and C, there is only one mode and the SSA is independent of m within the BR, being R = −1; see Appendix B for a rigorous proof.
In the macroscopic limit of interest in optics, we approximate q ν (l c + j c ) 1 and hence in Eq. ( 16) take e −iq ν l ≈ 1 and S 0 (l , m) ≈ S 0 (0, m), thence obtaining where T = l,l T (l, l ).For l c = 1, this reduces to a constant R = T s independent of m; if assuming j c = 0 in addition, T = −1/s, thereby recovering the result that R = −1.Equation (17) shows that R ν (m) can in general be written as a sum of exponentials, as proposed in Ref. [28].

C. Slab system
Suppose the slab consists of N actual layers of oscillators, viz.n = 1, . . ., N. It possesses two SRs: the left one that is comprised of the layers with 1 − l c n j c and the right one of layers with N + 1 − j c n N + l c .The remaining layers make up the BR; see Fig. 1.Note that the layers with n < 1 or n > N are virtual layers as described above.
The solution to the equation of motion, that is, Eq. ( 4), can be written as r(n) = (e/M ) m S(n, m)E (m), where S(n, m) is given by the Lippmann-Schwinger relation, Eq. ( 14), but with H and the corresponding scattering matrix T defined for the slab rather than the semi-finite system.
Since the slab has two SRs, let us write H = V l + V r , where V l and V r refer to the left SR (lSR) and the right SR (rSR), respectively.As such, V l (l, l ) = δ l,l (l ) for l, l j c and V r (l, l ) = δ l,l (l ) for l, l N + 1 − j c and all remaining elements of V l and V r vanish.One may also introduce T l/r = −V l/r (1 + S 0 V l/r ) −1 , which describes the scattering within the l/rSR.Clearly, T l/r (l, l ) vanishes if l or l or both lie in the r/lSR (and in the BR).Applying Eq. ( 14) and the expression of S 0 , Eq. ( 9), one obtains that for n in the BR Here m = N + 1 − m and n = N + 1 − n, which measure the distances to the rSR.Now we have three waves to each mode: the original waves excited by the electric field at layer m (i.e., e iq ν |n−m| ), the right-going waves (i.e., e iq ν (n+m) ), and the leftgoing waves (i.e., e iq ν (n+ m) ).The amplitudes R l/r ν (m) are given by Note that the scattering matrix T can be expressed in terms of S 0 , T l , and T r .Actually, by regrouping the terms in the series of T , one can show that T = Tl + Tr , where Tl/r split as Tl/r = T (1) l/r + T (2) l/r with T (1) l/r = T l/r (1 − S 0 T r/l S 0 T l/r ) −1 , T (2) l/r = T (1) l/r S 0 T r/l .(21) One may see that T (l ∈ l/rSR, l ) = Tl/r (l, l ) and Tl/r (l ∈ r/lSR, l ) = 0.As a result, which resembles Eq. ( 16).In the presence of dissipation, q ν carries a positive imaginary part and the elements of S 0 involved in Eq. ( 21) are negligible for sufficiently large N (see below).Then Tl/r reduces to that for a single SR and R l/r ν (m) reduces to the SSAs R ν (m).
Physically, Tl/r (l, l ) describes the process that an incident wave hits layer l and subsequently layer l ∈ lSR gets activated to generate right/left-going waves.According to which SR the layer l belongs, Tl/r (l, l ) naturally partitions into T (1) l/r (l, l ), which collects the contributions with l ∈ l/rSR, and T (2) l/r (l, l ), which collects the contributions with l ∈ r/lSR.If the incidence wave hits a layer in one SR, a layer in the other SR will get activated only if the SR being hit generates a wave that travels to the other layer and activates the latter.This explains the relation that T (2) l/r = T (1) l/r S 0 T r/l .Needless to say, T (1) l/r (l, l ) vanishes unless both l and l lie in the l/rSR, whereas T (1) l/r (l, l ) vanishes unless l lies in the l/rSR while l lies in the r/lSR.
It is worth noting that T (1/2) l/r (l, l ) can be broken down into an infinite series of subprocesses of waves bouncing back and forth between the SRs.This becomes clear if we extend Then T (1) l/r = T l/r + T l/r S 0 T r/l S 0 T l/r + • • • , where the first term represents scattering within one SR and therefore involves no bounces, but the additional terms represent scattering involving a number of bounces between the SRs.For example, the second term T l/r S 0 T r/l S 0 T l/r involves two bounces; it describes the process that an incidence wave is scattered (i.e., the factor T l/r ) by the l/rSR, which is then activated to generate waves that propagate (S 0 ) to the r/lSR, where such waves are scattered (T r/l ) and the SR is activated to generate waves that travel (S 0 ) to the l/rSR, where scattering (T l/r ) takes place again to generate further waves.Similar interpretation applies to the remaining terms and the terms contained in T (2) l/r .It is clear that the elements of S 0 involved in T (1/2) l/r describe the propagation of waves between the SRs; hence they scale as S ν e iq ν N , vanishing for N → ∞.
In the macroscopic limit that the thicknesses of the SRs are much thinner than any other length scales of the system (e.g., q −1 ν ), we may take S 0 (l, l ) ≈ S0 = ν S ν e iq ν N if l and l belong to different SRs.Then T (1) l/r = l,l T (1) l/r (l, l ) ≈ Tl/r (1 − S0 Tr/l S0 Tl/r ) −1 , where Tl/r = l,l T l/r (l, l ) and T (2) l/r = l,l If the slab is symmetric about its midplane, Tl = Tr = T and Together with S 0 (N, m) = S 0 (0, m), this leads to showing that R l ν and R r ν are given by the same function R ν (m) = ( T (1) S 0 (0, m) + T (2) S 0 (0, m))e −iq ν m (27) due to symmetry.Explicitly, For nearest-neighbor coupling, this bears a simple relation to the SSA given as R = T s [see Eq. ( 17)], This relation was also derived in Ref. [28] by a simple physical argument.While R is a constant independent of m, R does depend on m due to multiple reflections of waves between the SRs.It was shown that the effects due to multiple reflections are not captured by ABCs [28].More generally, R ν as given by Eq. ( 28) can be expressed in terms of the SSAs R ν as given by Eq. ( 17).Actually, Eq. ( 17) can be rewritten as R ν = μ M νμ S μ , where the matrix M is introduced with elements M νμ = T e i(q μ −q ν )m .Inverting this matrix gives S μ = ν (M −1 ) μν R ν , which can then be substituted into Eq.( 28) to obtain the desired relation.

III. OPTICAL RESPONSES
With the dielectric responses of the system described in the preceding section, we now examine its interaction with light.For simplicity, let us consider a light field E i (z) normally incident unto a slab [31].Then the discrete structure along the x-y plane can be replaced by a uniform distribution of dipoles.The electric field and the resulting polarization field are taken to be along the x direction.The polarization field is given by where we have restored the lattice constant a for clarity and δ(z) is the Dirac function.The electric field is given by [28] where k 0 = ω/c is the wave number of the incidence field [i.e., E i (z) = e ik 0 z ], with c being the speed of light in the vacuum.Writing E (z n ) as E (n), Eq. ( 31) produces This equation, combined with and introducing and a dimensionless parameter b = 2π ik 0 a e 2 /a which can be compactly written as (1 − bK )E = E i , where K, E, and E i are matrices with elements K (n, m), E (n), and E i (n), respectively.Obviously, E = (1 − bK ) −1 E i , from which the polarization field can be computed.
We have numerically solved Eq. ( 35) with S computed according to Eqs. ( 13) and (14), which are applied here to the slab system.An analytical solution for the case with l c = 1 and j c = 0 is reported in Appendix C. A few optical quantities, i.e., reflectance = |E (0) − 1| 2 , transmission = |E (N )| 2 , and absorption = 1 − reflectance − transmission, are obtained across a wide spectrum of frequencies.They are displayed in Fig. 2 for a thin slab (N = 50) and a thicker one (N = 300).In obtaining the results, nearest-neighbor coupling is assumed with L 1 /ω 2 0 = −0.5 so that ω 2 q /ω 2 0 = 1 − cos(qa) and S 0 is given by Eq. (10).Other parameters are as follows:  As seen in Fig. 2, for the thin slab the optical quantities exhibit strong oscillations against ω within the continuum ω < ω q (i.e., ω/ω 0 < √ 2); such oscillations are much less pronounced for the thick slab due to a much shorter oscillation period.The oscillations are due to the factor e iqN that signifies the bounces of waves between the SRs.The oscillation period is therefore given by ω = 2π Na/v g , where v g = ∂ω ∂q is the group velocity of the polarization waves.Note that 2π/ ω gives the time it takes for a polarization wave packet to traverse the slab.For the parameters given above, v g ≈ aω 0 / √ 2 for ω/ω 0 √ 2 and hence ω ≈ (π √ 2/N )ω 0 , which agrees nicely with what is seen in Fig. 2.
The other interesting feature to note in Fig. 2 is that, outside the continuum, where ω/ω 0 > √ 2, there exists a strong absorption peak.This is due to the formation of surface states, one located within each SR.These states locate the poles of the matrix T [and hence the SSA R ν by way of Eq. ( 15)].For thick slabs, the surface states at different SRs do not overlap and their frequencies are independent of slab thickness.There is no surface state for hard-wall surfaces ( j c = 0).For the parameters given in the above, the spatial profile of the surface states is shown in Fig. 3, where the dipole d (n) (upper panel) is plotted at the corresponding frequency together with the field E (lower panel).It is seen that d (n) is strongly concentrated within the SRs and quickly decays into the BR.With this, it follows from Eq. ( 32) that E (n) ≈ e ik 0 na + 2π ik 0 a −2 d (1)(e ik 0 na + e ik 0 (N−n)a ).Since Nk 0 a 1, E (n) is almost a constant as seen in Fig. 3.In the macroscopic limit, where the thickness of the SRs is negligible in comparison with other length scales of the system (i.e., q j c 1, k 0 j c 1, and N j c ), the fine details within the SRs cannot be probed.As such, one may extend the validity of Eq. ( 18), the expression of S(n ∈ BR) supposed to be valid only for n lying in the BR, to the entire slab including the SRs.Let us call the resulting quantity S SSA (n, m), which equals S(n, m) for n ∈ BR but differs for n ∈ SRs.In Fig. 4,  S(n, m) and S SSA (n, m) are plotted alongside for comparison.Since the difference only occurs within the SRs, one can approximate S(n, m) by S SSA (n, m) in the macroscopic limit, which underlies the SSA theory derived in Ref. [28].It is worth seeing that S SSA is also a general solution to Eq. ( 8).In Fig. 5, we display the optical quantities calculated using the exact S(n, m) alongside those by S SSA (n, m) for the thick slab (N = 300).They agree well with each other as long as ω lies far from the bound state.Close to the bound states, discrepancy occurs, which is expected: the bound states are localized within the SRs, where S SSA (n, m) differs from S(n, m).For the thin slab (N = 50), as shown in Fig. 6, good agreement occurs only for not so big ω (and hence q); as ω increases, a larger wave number is involved and the discrepancy becomes pronounced due to oscillations caused by multiple reflections which magnify the small mismatch between the S SSA and the exact SSA.As aforesaid, such oscillations are significant for thin slabs.Notwithstanding, the overall trend is still well captured by S SSA even at high frequencies.
Finally, the long-wavelength dielectric function of the infinite system can be obtained as ε(Q, ω) = 1 + 4πe 2 Ma 3 1 ω 2 Q − ω2 , with Q being the wave number.The dispersion of the polariton waves is then determined by ε(Q, ω) = (Q/k 0 ) 2 , which admits two pairs of roots, ±Q 1 and ±Q 2 , for each value of ω.Their expressions can be found in Ref. [28].By the conventional wisdom, this would beckon the ABCs in determining the optical responses of the system.However, as shown in the above and in Ref. [28], there is no need of ABCs at all.The optical properties are uniquely and fully determined with boundary effects encrypted in the SSA.

IV. SUMMARY
To summarize, we have studied a spatially dispersive model dielectric which by the conventional wisdom would require ABCs to uniquely determine its optical properties.Consistent with the SSA theory developed in Ref. [28], we have shown that no ABCs are needed and the macroscopic surface effects are completely captured by the SSAs, which we have calculated for the model.Universal analytical expressions are established for the SSAs and expressed in terms of the microscopic surface scattering matrix.Their various properties have been revealed.We showed that for nearest-neighbor coupling the SSAs do not depend on where the incident waves were excited [cf.Eq. ( 17)] but this is not so if coupling occurs at longer range.We also showed that, while they are defined for a semi-infinite system with a single surface, the same SSAs also determine the surface effects in the case of a slab, consistent with the theory derived in Ref. [28].
While the calculations are explicitly done for a model dielectric, we expect them to be widely applicable.This is so because, as is clear from our derivations, the main results rely on only two assumptions: (1) translational symmetry of the infinite system so that S 0 can be expressed as a sum of outgoing (quasi)plane waves and (2) the Lippmann-Schwinger relation, which is valid for any (linear) wave scattering theory and defines the SSA structure.As long as these assumptions apply, the main results, given in Eqs. ( 15), ( 17), (18), and (28), should remain valid.Indeed, these expressions do not involve any quantities exclusive to the model system.In addition, the formalism can be extended to boundaries of nonplanar geometries such as a cylinder or a sphere.
Finally, we should mention that the results can be experimentally studied as discussed in Ref. [28].For example, the SSA can be extracted from the variations of the optical reflection or transmission as a function of the thickness of a semiconducting slab, as implied in Eq. ( 29).This would require varying the thickness of the slab while leaving the surface morphology intact, which is routinely achieved using a growth rate gradient in moleular beam expitacy [32].Alternatively, one may study the optical transmission or reflection by two parallel slabs spaced by a dielectric, such as vacuum, gas, or liquid, having a thickness which can be tuned, for example, by a microelectromechanical system [33].
where we have left out the self-sustained part and in obtaining the second equality we have used the already proven result that λ = 0. We thus conclude that S(k, m) ≡ 0 for k < 1 − l c and m > 0, as claimed.Physically, Eq. (A5) says that the polarization on the virtual layers acts as a direct source for the polarization on layers to the left; as the polarization on the virtual layers is zero there is no polarization to the left.
The above result implies that a stack of at least l c virtual layers cuts a system into two decoupled subsystems.In the main text, we stated that R = −1 for the nearestneighbor coupling (l c = 1) with a hard-wall surface ( j c = 0).This corresponds to the well known fact that a wave changes its phase by π when reflected by a hard-wall surface, on which the wave has a node.Here we formally prove this.Let us start with Eq. ( 16), which for the present case reads R = T (0, 0)S 0 (0, m)e −iqm = T (0, 0)S 0 (0, 0).For the hardwall surface, T (0, 0)S 0 (0, 0) = −1, thus proving the claim.

FIG. 1 .
FIG.1.Slab divided into a left surface region (lSR), a right surface region (rSR), and a bulk region (BR).The size of blue disks indicates the value of , which vanishes in the BR.Gray dots represent virtual oscillators on which → ∞.

FIG. 2 .
FIG. 2. Reflectance, transmission, and absorption for slabs (solid black lines for N = 300 and circled pink lines for N = 50) with nearest-neighbor coupling.Parameters are given in the text.

FIG. 5 .
FIG. 5. Reflectance and transmission computed by S (exact) and those by S SSA (SSA) for N = 300.
APPENDIX B: SSA WITH l c = 1 AND j c = 0

APPENDIX C: OPTICAL PROPERTIES WITH l c = 1
AND j c = 0