Tunneling of multi-Weyl semimetals through a potential barrier under the influence

of magnetic fields

Ipsita Mandal
Faculty of Science and Technology, University of Stavanger, 4036 Stavanger, Norway and
Institute of Nuclear Physics, Polish Academy of Sciences, 31-342 Kraków, Poland

Sandip Bera
Indian Institute of Technology, Kharagpur, Kharagpur 721302, India

CONTENTS

I. Model and Formalism for 2D 1

II. Model and Formalism for 3D 3

III. Numerical results and discussions 4

A. 2d case 5
1. type one 5
B. new1:line number 593 in latex 6
1. type two 6
C. 3d case 7
D. Identifying the resonances with quasibound states 8
1. Quasibound states within a 2D static quantum well 8
2. Quasibound states within a 3D static quantum well 11
E. Relation between different components of shot noise 13

A. DETAILED CALCULATION OF FLOQUET SCATTERING MATRIX FOR 2D 14

B. DETAILED CALCULATION OF FLOQUET SCATTERING MATRIX FOR 3D 16

References 27

I. MODEL AND FORMALISM FOR 2D

The model Hamiltonian for the quadratic band crossing point in 2d system is described by the Hamiltonian [? ?
] H = H + V (r, t) where

1
2px py σx + (p2x − p2y )σz


H=− (1.1)
2µ

V (r, t) = V0 + V1 cos(ωt), − L/2 ≤ x ≤ L/2 (1.2)

=0 otherwise

where µ is the effective mass and V0 is the static potential. The time-dependent potential V (r, t) is spanned only
along x direction with driving amplitude V1 and frequency ω. So, the time-dependent Schrödinger equation can be
written
~2
σx ∂x ∂y + σz (∂y2 − ∂x2 ) ψ(x, y, t) + σ0 V (x, t)ψ(x, y, t)


i~∂t ψ(x, y, t) = − (1.3)
2µ
 2

The Floquet theorem asserts [1] that Eq.(1.3) has solutions of the form ψ(x, y, t) = e−iEn t/~ ψ(x, y, t) where
En is the Floquet eigen energy. The potential is translation invariant in the y-direction therefore a plane
wave component eiky y can be assumed with ky conserved during transmission .So, Floquet spinor ψ(x, y, t) =
eiky y (ψ1 (x, t), ψ2 (x, t))T . The mathematical form of the spinor inside and outside of the quantum well can be
written as

X
e−iEn t/~ Ain eiknx x + Aon e−iknx x γ1


ψ1 (x, t) = x < −L/2 (1.4)
n
 
X X V1
e−iEn t/~ γ2 am eiqm x + bm e−iqm x Jn−m


= Θ(E − V0 ) (1.5)
n m
~ω
 
iqm x V1
−iqm x


+ γ3 am e + bm e Jn−m Θ(V0 − E), − L/2 ≤ x ≤ L/2
~ω
X
e−iEn t/~ Bni e−iknx x + Bno eiknx x γ1


= x > L/2 (1.6)
n

X
ky
ψ2 (x, t) = e−iEn t/~ −Ain eiknx x + Aon e−iknx x γ1 x < −L/2 (1.7)
n
knx
 
X X V1 ky
e−iEn t/~ −am eiqm x + bm e−iqm x Jn−m


= γ2 Θ(E − V0 )+ (1.8)
n m
~ω qm
 
iqm x −iqm x

qm
V1
am e − bm e Jn−m γ3 Θ(V0 − E), − L/2 ≤ x ≤ L/2
ky
~ω
X
ky
= e−iEn t/~ Bni e−iknx x − Bno sn eiknx x γ1 x > L/2 (1.9)
n
knx

knx qm ky
where γ1 = √ , γ2 = √ and γ3 = √ are the normalization factors. Here, Ain and Aon are
ky2 +knx
2 ky2 +qm
2 ky2 +qm
2

the probability amplitudes of the incoming and outgoing waves from the left and Bni and Bno are those from right
respectively. am and bm are the amplitudes of wave function in the oscillating region that can be determined by the
q conditions. En =qEF + n~ω is the Floque eigenenergy of the nth order Floque state with n = 0, ±1 ± 2.. .
boundary
2µ|Em −V0 |
knx = 2µE 2
~2 − ky , qm =
n
~2 − ky2 are the wave vectors without and with oscillation potential respectively.
Jn (x) are nth-order first kind Bessel functions, which only exists in the oscillating region [2]. ky remains same in
all regions because of the nature of the potential.
The wave function must be continuous at the boundaries x = ±L/2. After doing some algebra, we obtain the
following matrix equation (see the Appendix A (A)):
 o  i 
An P Am
= S
m nm
Bon Bim
The matrix in above expression reflect a connection between nth Floquet channel to mth Floquet channel.
Here, this scattering matrix is known as Floquet scattering
q matrix. Considering the real current flux, the Floquet
Re(knx )
scattering matrix elements [3] are sαβ (En , Em ) = Re(k mx )
Snmαβ . Here s(En , Em ) is a matrix but sα,β (En , Em )
is scalar quantity. The general formula for Floquet scattering matrix becomes
q
s(En , Em ) = kkmx
nx
Snm
rnm t0nm
   
s (E , E ), sLR (En , Em )
= LL n m = 0
sRL (En , Em ), sRR (En , Em ) tnm rnm
with rnm and tnm are the reflection and transmission amplitudes from the mth Floquet channel to the nth Floquet
0
channel respectively. rnm and t0nm are the corresponding backward amplitudes. Total transmission and reflection
probability
T = n |t0n |2 = |sRL (EF , En )|2 R = n |r0n |2 = |sLL (EF , En )|2
P P

Zero-frequency nonadiabatic pumped shot noise at low temperature can be expressed as [4–8]
 3

∞ ∞ 2
e2 Mαβγδ (E, Em , En , Ep ) [f0 (En ) − f0 (Em )]
Z X X
Nαβ = ,
h 0 m,n,p=−∞
2
γ,δ=L,R

Mαβγδ (E, Em , En , Ep ) = sαγ (E, Em ) sαδ (E, Em ) s∗βδ (Ep , Em ) sβγ (Ep , En ).
∗
(1.10)

II. MODEL AND FORMALISM FOR 3D

Total solution in three different regions are

X
φ(x, y, x, t) = e−iEn t/~ eiky y eikz z ψL x < −L/2 (2.1)
n
 
V1
ψM Jn−m − L/2 ≤ x ≤ L/2 (2.2)
~ω
ψR x ≥ L/2 (2.3)

where

fil11 fil21 fol11 fol21

       
f  f  f  f 
ψL = Ain,1 eikn x  il12  + Ain,2 eikn x  il22  + Aon,1 e−ikn x  ol12  + Aon,2 e−ikn x  ol22  (2.4)
fil13 fil23 fol13 fol23
fil14 fil24 fol14 fol24
fim11 fim21 fom11 fom21
       
X f f f f
ψM = αm,1 eiqm x  im12  + αm,2 eiqm x  im22  + βm,1 e−iqm x  om12  + βm,2 e−iqm x  om22  (2.5)
   
fim13 fim23 fom13 fom23
m
fim14 fim24 fom14 fom24
fil11 fil21 fol11 fol21
       
f  f  f  f 
o
ψR = Bn,1 eikn x  il12  + Bn,2
o
eikn x  il22  + Bn,1
i
e−ikn x  ol12  + Bn,2
i
e−ikn x  ol22  (2.6)
fil13 fil23 fol13 fol23
fil14 fil24 fol14 fol24

(kn + iky )(kz + χn ) i(3kz + χn ) i kn2 + ky2 − 2kz2 − 2kz χn 1
fil11 = − ; fil12 = √ ; f il13 = − √ ; fil14 =
n1 (kn − iky )2 3n1 (kn − iky ) 3n1 (kn − iky )2 n1


(kn + iky )(χn − kz ) i(χn − 3kz ) i kn2 + ky2 − 2kz2 + 2kz χn 1
fil21 = ; fil22 = − √ ; fil23 = − √ ; fil24 = ;
n2 (kn − iky )2 3n2 (kn − iky ) 3n2 (kn − iky )2 n2

(kn − iky )(kz + χn ) i(3kz + χn ) i kn2 + ky2 − 2kz2 − 2kz χn 1
f0l11 = 2
; f0l12 = − √ ; f0l13 = − √ ; f0l14 =
n1 (kn + iky ) 3n1 (kn + iky ) 3n1 (kn + iky )2 n1
2 2 2


(kn − iky )(χn − kz ) i(χn − 3kz ) i kn + ky − 2kz + 2kz χn 1
f0l21 =− ; f0l22 = √ ; f0l23 = − √ ; f0l24 =
n2 (kn + iky )2 3n2 (kn + iky ) 3n2 (kn + iky )2 n 2

2


(qm + iky )(kz + χm ) i(3kz + χm ) i qm + ky2 − 2kz2 − 2kz χm 1
fim11 =− ; fim12 = √ ; fim13 = − √ ; fim14 =
n11 (qm − iky )2 3n11 (qm − iky ) 3n11 (qm − iky )2 n11
2


(qm + iky )(χm − kz ) i(χm − 3kz ) i qm + ky2 − 2kz2 + 2kz χm 1
fim21 = 2
; fim22 = − √ ; fim23 = − √ ; fim24 =
n22 (qm − iky ) 3n22 (qm − iky ) 3n22 (qm − iky )2 n22
 4

2


(qm − iky )(kz + χm ) i(3kz + χm ) i qm + ky2 − 2kz2 − 2kz χm 1
f0m11 = ; f 0m12 = − √ ; f0m13 = − √ ; f0m14 =
n11 (qm + iky )2 3n11 (qm + iky ) 3n11 (qm + iky )2 n11
2


(qm − iky )(χm − kz ) i(χm − 3kz ) i qm + ky2 − 2kz2 + 2kz χm 1
f0m21 =− ; f0m22 = √ ; f0m23 = − √ ; f0m24 =
n22 (qm + iky )2 3n22 (qm + iky ) 3n22 (qm + iky )2 n22

q p q p
2 2
3 (2χn ) χn (kz + χn ) 3 (2χn ) χn (χn − kz )
n1 = ; n2 = ; (2.7)
kn2 + ky2 kn2 + ky2
q p q p
2 2
3 (2χm ) χm (kz + χm ) 3 (2χm ) χm (χm − kz )
n11 = 2 + k2
; n22 = 2 + k2
;
qm y qm y
q q q q
2µEn 2µ|Em −V0 |
kn = ~2 − ky2 − kz2 , qm = ~2 − ky2 − kz2 , χn = kn2 + ky2 + kz2 ,χm = 2 + k2 + k2 ,
qm y z
The wave function must be continuous at the boundaries x = ±L/2. After doing some algebra, we obtain the
following matrix equation (see the Appendix B (B)):

 0   i 
A1n A1m
A02n  X Ai 
2m 
B 0  = Snm  i 
B1m
1n m
0 i
B2n B2m
q
Re(χn )
s(En , Em ) = Re(χm ) Snm

s11 (En , Em ), s12 (En , Em ) s13 (En , Em ) s14 (En , Em ) r11nm r12nm t011nm t012nm
   
s21 (En , Em ), s22 (En , Em ) s23 (En , Em ) s24 (En , Em ) r21nm r22nm t021nm 0
t22nm 
= = 0 0
s31 (En , Em ), s32 (En , Em ) s33 (En , Em ) s34 (En , Em )  t11nm t12nm r11nm r12nm 
0 0
s41 (En , Em ), s42 (En , Em ) s43 (En , Em ) s44 (En , Em ) t21nm t22nm r21nm r22nm

Total transmission and reflection probability

X 1
|t0n |2 = |s11 (EF , En )|2 + |s12 (EF , En )|2 + |s21 (EF , En )|2 + |s22 (EF , En )|2


T =
n
2

X 1
|r0n |2 = |s31 (EF , En )|2 + |s32 (EF , En )|2 + |s41 (EF , En )|2 + |s42 (EF , En )|2


R=
n
2

∞ ∞ 2
e2 Mαβγδ (E, Em , En , Ep ) [f0 (En ) − f0 (Em )]
Z X X
Nαβ = ,
h 0 m,n,p=−∞
2
γ,δ=L,R

Mαβγδ (E, Em , En , Ep ) = sαγ (E, Em ) sαδ (E, Em ) s∗βδ (Ep , Em ) sβγ (Ep , En ).
∗
(2.8)

III. NUMERICAL RESULTS AND DISCUSSIONS

In this section, we show the plots of transmission coefficients and shot noise for some representative parameter
values of the systems. We compare our results with other systems like graphene, and also interpret the physical
implications of the results that we obtain.
 5

1.0 1.0
0.00095
0.8 0.00096 0.8
0.00097

0.6 0.6
R

T
0.4 0.4
0.00095
0.00096
0.2 0.2 0.00097

0.0 0.0
4.6 4.7 4.8 4.9 5.0 5.1 5.2 5.3 4.6 4.7 4.8 4.9 5.0 5.1 5.2 5.3
Ef Ef
(a) (b)

FIG. 1. 2d case: Total Floquet reflection and transmission coefficient for different sets of ky (0.00095Å−1 (blue) ,0.00096Å−1
(green) and 0.00097Å−1 (red)), through time dependent quantum well. Sharp Fano resonance can see in the figures (1a) and
(1b) at a specific value of EF . The numerical parameters are ~ω = 4.0mev, L = 3000Å, V0 = −10 mev and V1 = 1 mev.

0.00095 1.2
Derivative of shot noise (10-2e2h)

0.40 0.00096
1.0
Shot noise (10-2e2ω)

0.00097
0.38
0.8
0.36
0.6
0.34

0.32 0.4
0.00095
0.30 0.2 0.00096
0.00097
0.28 0.0
4.80 4.85 4.90 4.95 5.00 5.05 4.80 4.85 4.90 4.95 5.00 5.05
Ef Ef
(a) (b)

FIG. 2. 2d case: (2a) Pumped shot noise, and (2b) differential pumped shot noise driven by an oscillating potential well as
a function of EF for different ky . An inflection occurs in the pumped shot noise corresponding to the Fano resonance in the
transmission. Sharp resonance can be seen at the inflection energy in panel

A. 2d case

1. type one

The numerical results of reflection and transmission coefficient as a function of the EF are shown in figures (1)
and (3). The summation limit depends on N, number of sidebands that need to be included is determined by the
strength of the oscillation from N > V1 /~ω. In our numerical calculation, we choose N = 2 which is satisfied
the relation V1 << ~ω. Fano resonance [1, 9, 10] occurs when the first order Floquet sideband overlaps with
the quasibound level, which shows clearly in the figures (1) and (3). Instead of symmetric [11, 12] , We have a
 6

1.0 1.0
0.00095
0.0008
0.8 0.8

0.6 0.6

T
0.4 0.4

0.2 0.00095 0.2

0.0008

0.0 0.0
4.0 4.5 5.0 5.5 6.0 4.0 4.5 5.0 5.5 6.0
Ef Ef

(a) (b)

FIG. 3. 2d case: Total Floquet reflection (fig. 3a) and transmission (fig. 3b) coefficient with EF for ky = 0.00095Å−1 (red)
and, ky = 0.0008Å−1 (green).

Derivative of shot noise (10-2e2h)

Derivative of shot noise (10-2e2h)

0.38 5.1
1.195 5.
Shot noise (10-2e2ω)

Shot noise (10-2e2ω)

0.36 3.4
1.19 4.3

0.34 1.7

1.185 3.7
0.32 -0.077
5.77 5.78 5.79 5.80 5.81 4.75 4.80 4.85 4.90
Ef Ef

(a) (b)

FIG. 4. 2d case: Pumped shot noise, and differential pumped shot noise driven by an oscillating potential well as a function
of EF for ky = 0.0008Å−1 (fig. 4a) and ky = 0.00095Å−1 (fig. 4b)

asymmetric Fano resonance[13].

Fano resonance (FN) is clearly visible for ky = 0.00095Å−1 , 0.00096Å−1 and 0.00097Å−1 from the reflection and
transmission curve in figure (1). The evidence of the Fano resonance will increase by the shot noise and differential
shot noise, which are shown in figure (II). Figure (1) shows the Fano resonance point (FNP) increase with ky , which
is not consistent with the graphene results, as in graphene FNP decrease with increase of ky . The bandwidth of the
resonance curves are almost constant with ky . We check this feature for higher values of ky (up to 0.0015Å−1 ). Sim-
ilar to reflection and transmission amplitude, the value of the shot noise and differential shot noise are inconsistent
with graphene. Figure (2b) suggest, radius of curvature of shot noise is positive. But, it is negative for graphene.
The FN is weaker in lower limit of ky which is consistent with graphene.

B. new1:line number 593 in latex

The value of the shot noise and differential shot noise are inconsistent with graphene. Figure (2b) suggest,
the derivative of the differential shot noise changes sign from negative to positive when passing through the Fano
resonance point. Hence, curvature of the shot noise goes negative to positive if we go from left to right. So, Fano
resonance point on shot noise curve is an type of inflection point. Whereas for graphene, curvature of the shot noise
changes sign from positive to negative.

1. type two

Figure (3) shows the Fano resonance point (the value of EF at which transmission or reflection curve has a sharp
peak) decrease with increase of ky , which is consistent with the graphene. But, the value and sequence of the
 7

1.0 1.0 1.0 1.0 1.0 1.0

0.8 0.8 0.8 0.8 0.8 0.8

0.6 0.6 0.6 0.6 0.6 0.6

T

T
0.4 0.4 0.4 0.4 0.4 0.4

0.2 0.2 0.2 0.2 0.2 0.2

0.0 0.0 0.0 0.0 0.0 0.0

3 4 5 6 7 3 4 5 6 7 3 4 5 6 7 3 4 5 6 7 4 5 6 7 8 3.5 4.0 4.5 5.0 5.5 6.0
Ef Ef Ef Ef Ef Ef

(a) (b) (c) (d) (e) (f)

FIG. 5. 2d case: Figures (5a),(5a),(5b),(5c),(5d),(5e) and (5f) shows the total Floquet transmission coefficient at v0 = −10
for ky (0.00077Å−1 ,0.000978Å−1 , 0.00086Å−1 ,0.00087Å−1 ,0.00092Å−1 and, 0.00095Å−1 ) respectively.

1.0 1.0 1.0 1.0

0.8 0.8 0.8 0.8

0.6 0.6 0.6 0.6

T

T
0.4 0.4 0.4 0.4

0.2 0.2 0.2 0.2

0.0 0.0 0.0 0.0

4 5 6 7 8 9 4 5 6 7 8 9 4 5 6 7 8 9 4 5 6 7 8 9
Ef Ef Ef Ef
(a) (b) (c) (d)

FIG. 6. 2d case: Figures (6a),(6b),(6c), and (6d) shows the total Floquet transmission amplitude at ky = 0.001Å−1 for
v0 = −5, v0 = −10, v0 = −15 and v0 = −20 mev respectively.

bandwidth with ky are not consistent with the graphene. The resonance peak for ky = 0.0008 is very sharp(small
bandwidth ) than ky = 0.00095. Nature of the differential shot noise in figure (4a) is completely different to the
previous one even with graphene. Here, derivative has a sharp pick at FNP. Extremely narrow bandwidth Fano-type
resonance can be useful for advanced dynamic reconfigurable devices [14].
So, based on bandwidth and the nature of the curves we can conclude, there are two different type of Fano
resonances. Region for different types of Fano resonance are shown in figure (5). In first case Fano resonance
gradually distorts as the ky decrease from 0.00095Å−1 (type one). At the same time small bandwidth type Fano
resonance will appear in the system (type two). Finally, both types of Fano resonance vanish in the limit ky ≤
0.00077Å−1
In figure (6) shows how Fano resonance depends on the static potential. The Fano resonance is completely
destroyed at v0 = −15 mev. But, Fano resonance rebuilds in the system at v0 = −20 mev. This whole process will
repeat with v0 . To check, we extent our calculation up to v0 = −45 mev. So, from result we can conclude that,
Fano resonance has a periodicity with unequal interval. Transmission resonances in figure (6a) at v0 = −5 mev
show two peak which means there are two static bound states [9].

C. 3d case

Similar to 2D, the Fano resonance point increase with ky for fixed value of kz and vice versa. To find the limit,
we extends our calculation over a large interval ( from 0.0002Å−1 to 0.0015Å−1 ) of ky and kz . But, we did not find
second type of Fano resonance. However, first type of FN exists over this interval, which is shown in figure (9).
This feature is completely different with 2D even with any other systems. Figure (8b) shows radius of curvature of
the shot noise is negative, which is consistent with graphene. Figure (11) shows how Fano resonance depends on
the static potential. This feature is similar with 2D.
 8

1.0 1.0
0.00095 0.00095
0.00096 0.00096
0.8 0.8
0.00097 0.00097

0.6 0.6
R

R
0.4 0.4

0.2 0.2

0.0 0.0
8.0 8.2 8.4 8.6 8.8 8.0 8.2 8.4 8.6 8.8
Ef Ef
(a) (b)

FIG. 7. 3d case: Total Floquet reflection and transmission coefficient in 3D system for different sets of ky (0.00095Å−1
,0.00096Å−1 and 0.00097Å−1 ) with kz = 0.00095Å−1 , through time dependent quantum well. Sharp Fano resonance can
see in the figures (7a) and (7b) at a specific value of EF .

0.70 12
0.00095 0.00095
Derivative of shot noise (10-2e2h)

0.00096 0.00096
10
Shot noise (10-2e2ω)

0.65 0.00097 0.00097

8
0.60
6

0.55 4

2
8.2 8.3 8.4 8.5 8.6 8.2 8.3 8.4 8.5 8.6
Ef Ef
(a) (b)

FIG. 8. 3d case: (8a) Pumped shot noise, and (8b) differential pumped shot noise driven by an oscillating potential well as
a function of EF for different ky at kz = 0.00095Å−1 .

D. Identifying the resonances with quasibound states

1. Quasibound states within a 2D static quantum well

We first consider the quasibound states within the static quantum well. The two-component fermion wavefunction
T
for a 2D static quantum well ψ(x, y) = eiky y (ψ1 (x), ψ2 (x)) can be decomposed as:
 9

1.0
1.0

0.8
0.8

0.6 0.6
T

T
0.4 0.4

0.00075
0.00075
0.2 0.2 0.00085
0.00085
0.00090
0.00090
0.0 0.0
7.0 7.5 8.0 8.5 5.5 6.0 6.5 7.0
Ef Ef
(a) (b)

FIG. 9. 3d case: Total Floquet reflection and transmission coefficient in 3D system for different sets of ky (0.00095Å−1
,0.00085Å−1 and 0.00075Å−1 ) with kz = 0.00095Å−1 and kz = 0.0007Å−1 respectively. Sharp Fano resonance can see in the
figures (9a) and (9b) at a specific value of EF .

1.0 0.00075
Derivative of shot noise (10-2e2h)

12
0.00085
0.0009
0.8
Shot noise (10-2e2ω)

10

0.6 8

0.4 6
0.00075
0.00085
4
0.2 0.00090

5.4 5.6 5.8 6.0 6.2 6.4 6.6 6.8 5.6 5.8 6.0 6.2 6.4 6.6 6.8
Ef Ef
(a) (b)

FIG. 10. 3d case: (10a) Pumped shot noise, and (10b) differential pumped shot noise for different sets of ky (0.00095Å−1
,0.00085Å−1 and 0.00075Å−1 ) at kz = 0.00075Å−1 .

kx x
for x < − L2

γ
 h 1 A1 e i
γ2 A2 eiqx + B2 e−iqx Θ(E + V0 ) + γ3 A2 eiqx + B2 e−iqx Θ(−E − V0 ) L L



ψ1 (x, t) = for − 2 ≤x≤ 2 ,

γ1 A3 e−kx x L

for x > 2
(3.1)
 10

1.0 1.0 1.0 1.0

0.8 0.8 0.8 0.8

0.6 0.6 0.6 0.6

T

T
0.4 0.4 0.4 0.4

0.2 0.2 0.2 0.2

0.0 0.0 0.0 0.0

8 9 10 11 12 13 8 9 10 11 12 13 8 9 10 11 12 13 8 9 10 11 12 13
Ef Ef Ef Ef
(a) (b) (c) (d)

FIG. 11. 3d case: Figures shows the total Floquet transmission amplitude at ky = kz = 0.001Å−1 for v0 = −5, v0 = −10,
v0 = −15 and v0 = −20 mev respectively.


γ1 ky
hikx A1 e


kx x
i for x < − L2
γ2 ky γ3 q
−A2 eiqx + B2 e−iqx Θ(E + V0 ) + −iqx L L



ψ2 (x, t) = q ky A2 eiqx − B2 e Θ(−E − V0 ) for − 2 ≤x≤ 2 ,

− γ1 ky A e−kx x

for x > L
ikx 3 2
(3.2)
where
r r
1 1 1 2µE 2 µ |E + V0 |
γ1 = q , γ2 = q , γ3 = q , kx = − 2 + ky2 , q= − ky2 , (3.3)
1+
ky2 ky2
1+ q2 ~ ~2
kx2 1+ q2 ky2

and −V0 is the static potential inside the well. By continuity of the wave function at the boundaries of the
quantum well, the secular equation in the limit E > −V0 for the eigenenergy can be written as

kx L qL qL
γ1 e− 2 −γ2 e− 2 −γ2 e 2 0
qL qL kx L
0 γ2 e 2 γ2 e− 2 −γ1 e− 2
kx L qL qL
ζ= γ1 e− 2 ky γ2 e− 2 ky γ2 e 2 ky =0 (3.4)
ikx q − q 0
qL qL kx L
γ2 e 2 ky γ2 e − 2 k y γ1 e− 2 ky
0 − q q ikx

~2 ky2 ~2 ky2
Roots of E for this equation are the quasibound (Eb ) state energies. 2µ − V0 and 2µ are the lower and upper
~2 ky2
limit of E respectively. is a positive quantity. Numerical results of the secular equation are given in the small
2µ
black dots in fig. (12). Shot noise and differential shot noise as a function of the Fermi energy are shown in Figure
(2). The differential shot noise curve (blue) for ky = 0.00095Å−1 have a sharp dip at Ef = 4.844meV . Now from
the relation EF ano − ~ω = Eb , Eb is 0.844 meV (~ω = 4 meV). At the same condition, the energy of the quasibound
state (from numerical analysis in fig. (12)) is Eb = 0.844 meV. Hence, quasibound state exists inside the quantum
well. Figure (12) suggests, there are many quasi bound states for a given ky . The sideband energy gap ~ω will
increase if we increase the driving frequency. In this situation, higher quasibound states can be activated into the
transport process. The above mechanism suggests, particle in the propagating state can emit photons with energy
~ω and drop into the quasibound state. The dependency of EF ano on Eb in the above formula is consistent with
graphene. So the Fano resonance in presence of oscillating potential is completely controlled by the energy of the
static quasibound state.
Not for paper
q
• In [9]: Here the quantum well is a 1D system. So, outside the quantum well, kx = 2µE ~1 (for statics case
) which is always real. therefore E > 0, solution e±kx x and, the state inside the quantum well is known as
bound state.
 11

• However,
q in our case (2D) or in ref. [5, 12], the system quantum well is a 2D. So, outside the quantum well,
kx = − 2 ~µ2E + ky2 (for statics case )

q
1. Now if 2 ~µ2E > ky2 . kx is imaginary. Then solution is e±iχx . Where χ = 2µE
~2 − ky2 . The state inside
the quantum well is known as bound state.
2. There is another limit, 2 ~µ2E < ky2 . kx is real. Then solution is e±kx x . The state inside the quantum well
is known as quasi-bound state. Our analysis and in ref. [5, 12] are belong to this limit.

15

10

Eb / EFano-ℏω (meV)
5

-5

- 10
0.5 1.0 1.5 2.0
ky (10-3 Å-1)

FIG. 12. Black points are dispersion of the quasibound states Eb of the static quantum well with depth V0 and width L in
2d QBT. The red and blue points represents the value of EFano − ~ ω for type one and type two Fano resonance respectively.
The numerical parameters used for the driven well are: ~ ω = 4 meV, L = 3000 Å, V0 = 10 meV, and V1 = 1 meV.

2. Quasibound states within a 3D static quantum well

The two-component fermion wavefunction for a 3D static quantum well ψ(x, y) = eiky y+ikz z φ(x) can be decom-
posed as:

    

 fil11 fil21
1 kx x fil12  1 kx x fil22 

    
for x < − L2


 A 1 e   + A 2 e  



 fil13  fil23 

 f f
 il14  il24 


    
fim11 fim21 fom11 fom21 #



"

2 iqx fim12  2 iqx fim22  2 −iqx fom12  2 −iqx fom22 

        

 A 1 e   + A 2 e   + B1 e   + B 2 e   Θ(E + V0 )
fim13  fim23  fom13  fom23 





 f f f f
φ(x) = im14  im24  om14  om24  ,

 " g im11 g im21 gom11 gom21 #

 g  g  g  g 
  im12   im22   om12   om22 
+ A21 eiqx   + A22 eiqx   + B12 e−iqx   + B22 e−iqx   Θ(−E − V0 ) for − L
≤x≤ L


2 2



 gim13  gim23  gom13  gom23 

g g gom14 gom24


 im14  im24


 



 f ol11 f ol21

 f  f 
 3 −k x ol12 3 −k x ol22  L
B1 e x   + B2 e x  for x >

   
 2
f f

ol13 ol23 

   


fol14 fol24
(3.5)
 12

where


(−ikx + iky )(kz + χ1 ) i(3kz + χ1 ) i −kx2 + ky2 − 2kz2 − 2kz χ1 1
fil11 =− , fil12 =√ , fil13 =− √ , fil14 = ,
n1 (−ikx − iky )2 3n1 (−ikx − iky ) 3n1 (−ikx − iky )2 n1


(−ikx + iky )(χ1 − kz ) i(χ1 − 3kz ) i −kx2 + ky2 − 2kz2 + 2kz χ1 1
fil21 = , fil22 = − √ , fil23 =− √ , fil24 = ,
n2 (−ikx − iky )2 3n2 (−ikx − iky ) 3n2 (−ikx − iky )2 n2
(3.6)

(−ikx − iky )(kz + χ1 ) i(3kz + χ1 ) i −kx2 + ky2 − 2kz2 − 2kz χ1 1
fol11 = , fol12 = −√ , fol13 =− √ , fol14 = ,
n1 (−ikx + iky )2 3n1 (−ikx + iky ) 3n1 (−ikx + iky )2 n1


(−ikx − iky )(χ1 − kz ) i(χ1 − 3kz ) i −kx2 + ky2 − 2kz2 + 2kz χ1 1
fol21 = − , fol22 = √ , fol23 =− √ , fol24 =
n2 (−ikx + iky )2 3n2 (−ikx + iky ) 3n2 (−ikx + iky )2 n2
(3.7)

(q + iky )(kz + χ2 ) i(3kz + χ2 ) i q 2 + ky2 − 2kz2 − 2kz χ2 1
fim11 =− , fim12 =√ , fim13 =− √ , fim14 = ,
n11 (q − iky )2 3n11 (q − iky ) 3n11 (q − iky )2 m1


(q + iky )(χ2 − kz ) i(χ2 − 3kz ) i q 2 + ky2 − 2kz2 + 2kz χ2 1
fim21 = , fim22 = − √ , fim23 =− √ , fim24 = ,
n22 (q − iky )2 3n22 (q − iky ) 3n22 (q − iky )2 m2
(3.8)

(q − iky )(kz + χ2 ) i(3kz + χ2 ) i q 2 + ky2 − 2kz2 − 2kz χ2 1
fom11 = , fom12 = − √ , fom13 = − √ , fom14 = ,
n11 (q + iky )2 3n11 (q + iky ) 3m1 (q + iky )2 m1


(q − iky )(χ2 − kz ) i(χ2 − 3kz ) i q 2 + ky2 − 2kz2 + 2kz χ2 1
fom21 = − , fom22 = √ , fom23 =− √ , fom24 = ,
n22 (q + iky )2 3n22 (q + iky ) 3m2 (q + iky )2 m2
(3.9)

q q q q
2 23 (3kx2 + χ21 ) kx ky + ky2 + kz (kz + χ1 ) 2 23 (3kx2 + χ21 ) kx ky + ky2 + kz (kz − χ1 )



n1 = 2 , n2 = ,
(kx2 + ky ) (kx + ky )2
q p q p
2
2 3 χ2 χ2 (kz + χ2 ) 2 23 χ2 χ2 (χ2 − kz )
m1 = , m2 = ,
q 2 + ky2 q 2 + ky2
r r
2µE 2 µ |E − V0 | 2µE 2 µ |E − V0 |
kx = − 2 + ky2 + kz2 , q= − ky2 − kz2 , χ1 = , χ2 = . (3.10)
~ ~2 ~2 ~2
We define the g’s in a similar manner as the f ’s. We will not need the explicit expressions for the g’s because in
our numerics we consider E > −V0 .

kx l kx l −iql −iql iql iql

fil11 e− 2 fil21 e− 2 −fim11 e 2 −fim21 e 2 −f0m11 e 2 −f0m21 e 2 0 0
kx l kx l −iql −iql iql iql
fil12 e− 2 fil22 e− 2 −fim12 e 2 −fim22 e 2 −f0m12 e 2 −f0m22 e 2 0 0
kx l kx l −iql −iql iql iql
fil13 e− 2 fil23 e− 2 −fim13 e 2 −fim23 e 2 −f0m13 e 2 −f0m23 e 2 0 0
kx l kx l −iql −iql iql iql
fil14 e− 2 fil24 e− 2 −fim14 e 2 −fim24 e 2 −f0m14 e 2 −f0m24 e 2 0 0
ζ= iql iql −iql −iql kx l kx l =0
0 0 fim11 e 2 fim21 e 2 f0m11 e 2 f0m21 e 2 −f0l11 e− 2 −f0l21 e− 2
iql iql −iql −iql kx l kx l
0 0 fim12 e 2 fim22 e 2 f0m12 e 2 f0m22 e 2 −f0l12 e− 2 −f0l22 e− 2
iql iql −iql −iql kx l kx l
0 0 fim13 e 2 fim23 e 2 f0m13 e 2 f0m23 e 2 −f0l13 e− 2 −f0l23 e− 2
iql iql −iql −iql kx l kx l
0 0 fim14 e 2 fim24 e 2 f0m14 e 2 f0m24 e 2 −f0l14 e− 2 −f0l24 e− 2
(3.11)
 13

Eb / EFano-ℏω (meV)
0

-5

- 10
0.2 0.4 0.6 0.8 1.0 1.2
ky (10-3 Å-1)

(a) (b)

FIG. 13. Blue points are dispersion energy of the quasibound states Eb of the static quantum well with depth V0 and width
L in 3d QBT. The red points represents the value of EFano − ~ ω. Fig. (13b) shows the zoomed version of the quasibound
states at kz = 0.00096Å−1 . The numerical parameters used for the driven well are: ~ ω = 4 meV, L = 3000 Å, V0 = 10 meV,
and V1 = 1 meV.

Roots of E for this equation are the quasibound (Eb ) state energies. Numerical results of the secular equation
are given in the small black dots in fig. (13). Shot noise and differential shot noise as a function of the Fermi
energy are shown in Figure (II). The differential shot noise curve for ky = 0.00095Å−1 and kz = 0.00096Å−1 have
a sharp dip at Ef = 8.401 meV. Now from the relation EF ano − ~ω = Eb , Eb is 4.409 meV (~ω = 4 meV). At
the same condition, the energy of the quasibound state (from numerical analysis in fig. (13)) is Eb = 4.41 meV.
Hence, quasibound state exists inside the quantum well. Figure (13) suggests, there are many quasi bound states
for a given ky . The sideband energy gap ~ω will increase if we increase the driving frequency. In this situation,
higher quasibound states can be activated into the transport process. The above mechanism suggests, particle in
the propagating state can emit photons with energy ~ω and drop into the quasibound state. The dependency of
EF ano on Eb in the above formula is consistent with graphene. So the Fano resonance in presence of oscillating
potential is completely controlled by the energy of the static quasibound state.

E. Relation between different components of shot noise

• For 2D: N11 = N22 = −N12 = −N21

For 3D: I1L + I2L = −(I1R + I2R ) = 0 [4]

4
X 4
X
s∗1γ s1δ + s∗2γ s2δ = − s∗3γ s3δ + s∗4γ s4δ



(3.12)
γ,δ=1 γ,δ=1

4
X 4
X 4
X 4
X
N11 = g s∗1γ s1δ s∗1δ s1γ N12 = g s∗1γ s1δ s∗2δ s2γ N21 = g s∗2γ s2δ s∗1δ s1γ N22 = g s∗2γ s2δ s∗2δ s2γ
γ,δ=1 γ,δ=1 γ,δ=1 γ,δ=1
4
X 4
X 4
X 4
X
N13 = g s∗1γ s1δ s∗3δ s3γ N14 = g s∗1γ s1δ s∗4δ s4γ N31 = g s∗3γ s3δ s∗1δ s1γ N41 = g s∗4γ s4δ s∗1δ s1γ
γ,δ=1 γ,δ=1 γ,δ=1 γ,δ=1
4
X 4
X 4
X 4
X
N23 = g s∗2γ s2δ s∗3δ s3γ N24 = g s∗2γ s2δ s∗4δ s4γ N32 = g s∗3γ s3δ s∗2δ s2γ N42 = g s∗4γ s4δ s∗2δ s2γ
γ,δ=1 γ,δ=1 γ,δ=1 γ,δ=1
4
X 4
X 4
X 4
X
N33 = g s∗3γ s3δ s∗3δ s3γ N34 = g s∗3γ s3δ s∗4δ s4γ N43 = g s∗4γ s4δ s∗3δ s3γ N44 = g s∗4γ s4δ s∗4δ s4γ
γ,δ=1 γ,δ=1 γ,δ=1 γ,δ=1
 14

∞
e2
R∞ P [fo (En )−fo (Em )]2
where g = h o
dE 2
m,n,p=−∞
So,
4
X
s∗1γ s1δ + s∗2γ s2δ (s∗1δ s1γ + s∗2δ s2γ )


N11 + N12 + N21 + N22 = g
γ,δ=1

4
X
s∗3γ s3δ + s∗4γ s4δ (s∗3δ s3γ + s∗4δ s4γ )


N33 + N34 + N43 + N44 = g
γ,δ=1
4
X
s∗1γ s1δ + s∗2γ s2δ (s∗1δ s1γ + s∗2δ s2γ )


=g
γ,δ=1

4
X
s∗1γ s1δ + s∗2γ s2δ (s∗3δ s3γ + s∗4δ s4γ )


N13 + N14 + N23 + N24 = g
γ,δ=1
4
X
s∗1γ s1δ + s∗2γ s2δ (s∗1δ s1γ + s∗2δ s2γ )


= −g
γ,δ=1

4
X
s∗3γ s3δ + s∗4γ s4δ (s∗1δ s1γ + s∗2δ s2γ )


N31 + N41 + N32 + N42 = g
γ,δ=1
4
X
s∗1γ s1δ + s∗2γ s2δ (s∗1δ s1γ + s∗2δ s2γ )


= −g
γ,δ=1

So, (N11 + N12 + N21 + N22 ) = (N33 + N34 + N43 + N44 ) = − (N13 + N14 + N23 + N24 ) = − (N31 + N41 + N32 + N42 )

Appendix A: DETAILED CALCULATION OF FLOQUET SCATTERING MATRIX FOR 2D

The wave function,ψ must be continuous at the boundaries x = ±L/2. At x = −L/2

 
  X   V1
γ1 Ain e−iknx L/2 + Aon eiknx L/2 = γ2 am e−iqm L/2 + bm eiqm L/2 Jn−m Θ(E − V0 ) (A1)
m
~ω
 
X 
−iqm L/2 iqm L/2
 V1
γ2 am e + bm e Jn−m Θ(V0 − E)
m
~ω

 k  

y
X  V1 ky
−Ain e−iknx L/2
+ Aon eiknx L/2
γ1 = −am e −iqm L/2
+ bm e iqm L/2
Jn−m Θ(E − V0 ) γ2 (A2)
knx m
~ω qm
 
X  V1 qm
+ am e−iqm L/2 − bm eiqm L/2 Jn−m Θ(V0 − E) γ3
m
~ω ky

similarly at x = L/2

 
  X  V1
Bni e−iknx L/2 + Bno eiknx L/2 γ1 = am e iqm L/2
+ bm e −iqm L/2
Jn−m Θ(E − V0 )γ2 + (A3)
m
~ω
 
X
iqm L/2 −iqm L/2
 V1
am e + bm e Jn−m Θ(V0 − E)γ3
m
~ω
 15

 k  

y
X  V 1 ky
Bni e−iknx L/2
− Bno eiknx L/2
γ1 = −am e iqm L/2
+ bm e −iqm L/2
Jn−m γ2 (A4)
knx m
~ω qm
 
X
iqm L/2 −iqm L/2
 V1 qm
+ am e − bm e Jn−m γ3
m
~ω ky

ky
2(Ain ± Bni )e−iknx L/2 γ1 = (A5)
knx
   
X 1 1 −iqm L/2 1 1 iqm L/2 V1
γ2 ky ( + )e ± ky (− + )e (am ± bm )Jn−m Θ(E − V0 )
m
q m knx qm k nx ~ω
   
ky qm −iqm L/2 ky qm iqm L/2 V1
+ γ3 ( − )e ±( + )e (am ± bm )Jn−m Θ(V0 − E)
knx ky knx ky ~ω

For simplicity we define some matrices whose matrix elements are:

   
± knx −iqm L/2 knx iqm L/2 V 1 γ2
(M1s )nm = (1 + )e ± (1 − )e Jn−m (A6)
qm qm ~ω γ1

   
± knx qm −iqm L/2 knx qm iqm L/2 V1 γ3
(M2s )nm = (1 − )e ± (1 + )e Jn−m (A7)
ky2 ky2 ~ω γ1

(Mr )nm = 2e−iknx L/2 δn,m (A8)

 
± γ2 V1
(M1c )nm = e−i(knx ±qm )L/2 Jn−m (A9)
γ1 ~ω

 
γ3 −i(knx ±qm )L/2 V1
(M2c ± )nm = e Jn−m (A10)
γ1 ~ω

(Mi )nm = e−iknx L δn,m (A11)

C±
m = am ± bm (A12)

So from the above equation

Mr · (Ai ± Bi ) = M1s ± · C± Θ(E − V0 ) + M2s ± · C± Θ(V0 − E) (A13)

 Mr · Ai
+ (M− + M− − (M−
 + −1 −1
 + −1 −1
Ao = M+ − Mi · Ai
 
jc · (Mjs ) js ) jc · (Mjs ) js ) (A14)
2
 + −1 − −1 −
 + −1 − −1
 Mr · Bi
+ M+
 
jc · (Mjs ) − (Mjs ) + Mjc · (Mjs ) + (Mjs )
2
i i
≡ A MAA + B MAB (A15)
 16

 Mr · Ai
Bo = M− − −1 − −1
  + −1  +
 + −1
jc · (Mjs ) + (Mjs ) + Mjc · (Mjs ) − (Mjs ) (A16)
2
 +  + −1 − −1 − + −1 − −1
 Mr · Bi
− Mi · Bi
 
+ Mjc · (Mjs ) + (Mjs ) + Mjc · (Mjs ) − (Mjs )
2
≡ Ai MBA + Bi MBB (A17)

where
  j =1 for E > V0 and
 j i=
 2 for E < V0
Ao MAA MAB A
=
Bo MBA MBB Bi
 i
A
=S
Bi
For
 onth band the above
 i matrix equation becomes
An P A
= m Snm
Bon Bi

Appendix B: DETAILED CALCULATION OF FLOQUET SCATTERING MATRIX FOR 3D

At x = −L/2

fil11 fil21 fol11 fol21

       
f  f  f  f 
Ain,1 e−ikn L/2  il12  + Ain,2 e−ikn L/2  il22  + Aon,1 eikn L/2  ol12  + Aon,2 eikn L/2  ol22  (B1)
fil13 fil23 fol13 fol23
fil14 fil24 fol14 fol24

fim11 fim21 fom11 fom21

       
X f f f f
= αm,1 e−iqm L/2  im12  + αm,2 e−iqm L/2  im22  + βm,1 eiqm L/2  om12  + βm,2 eiqm L/2  om22 
   
fim13 fim23 fom13 fom23
m
fim14 fim24 fom14 fom24
 
V1
× Jn−m
~ω

At x = L/2

fim11 fim21 fom11 fom21

       
 
X f f f f V1
αm,1 eiqm L/2  im12  +αm,2 eiqm L/2  im22  +βm,1 e−iqm L/2  om12  +βm,2 e−iqm L/2  om22  ×Jn−m
   
fim13 fim23 fom13 fom23 ~ω
m
fim14 fim24 fom14 fom24

fil11 fil21 fol11 fol21

       
f  f  f  f 
o
= Bn,1 eikn L/2  il12  + Bn,2
o
eikn L/2  il22  + Bn,1
i
e−ikn L/2  ol12  + Bn,2
i
e−ikn L/2  ol22  (B2)
fil13 fil23 fol13 fol23
fil14 fil24 fol14 fol24

From equation (B1)

 ikn L ikn L 1 1

e 2 f0l11 e 2 f0l21 e− 2 ikn L fil11 e− 2 ikn L fil21 Aon,1

 ik2n L ikn L 1 1   Ao 
e f e 2 f0l22 e− 2 ikn L fil12 e− 2 ikn L fil22
 ikn L 0l12
  n,2 
ikn L 1 1  i  =
 e 2 f0l13 e 2 f0l23 e− 2 ikn L fil13 e− 2 ikn L fil23  An,1
ikn L
e 2 f0l14
ikn L
e 2 f0l24
1
e− 2 ikn L fil14
1
e− 2 ikn L fil24 Ain,1
 17

iqm L iqm L
 1 1

e− 2 iqm L fim11 αm,1 + e− 2 iqm L fim21 αm,2 + e 2 f0m11 βm,1 + e 2 f0m21 βm,2
iqm L iqm L
X  e− 12 iqm L f

α
1
+ e− 2 iqm L fim22 αm,2 + e 2 f0m12 βm,1 + e 2 f0m22 βm,2
  
V1
 − 1 iqm L im12 m,1  × Jn−m (B3)

1 iqm L iqm L

m
e 2 fim13 αm,1 + e− 2 iqm L fim23 αm,2 + e 2 f0m13 βm,1 + e 2 f0m23 βm,2  ~ω
1 1 iqm L iqm L
e− 2 iqm L fim14 αm,1 + e− 2 iqm L fim24 αm,2 + e 2 f0m14 βm,1 + e 2 f0m24 βm,2

 
X V1
A0n,1 = Jn−m × (B4)
m
~ω
p 1
(kn − qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz − χn ) + iky (χn + χm ))
αm,1 p
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )

p 1
(kn − qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (kz − χm ) + qm (kz − χn ) + iky (χn − χm ))
+αm,2 p (B5)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )

p 1
(kn + qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz + χm ) + qm (χn − kz ) + iky (χn + χm ))
+βm,1 p (B6)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )

p 1
(kn + qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz − χm ) + qm (χn − kz ) + iky (χn − χm ))
+βm,2 p (B7)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
⇒

X
A0n,1 = nm
αm,1 v11 nm
+ αm,2 v12 nm
+ βm,1 v13 nm
+ βm,2 v14 (B8)
m

 
X V1
A0n,2 = Jn−m × (B9)
m
~ω
p 1
(kn − qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
− αm,1 p
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )

p 1
(kn − qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz − χm ) + qm (kz + χn ) − iky (χn + χm ))
−αm,2 p (B10)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )

p 1
(kn + qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (−kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
+βm,1 p (B11)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )

p 1
(kn + qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (kn (χm − kz ) + qm (kz + χn ) + iky (χn + χm ))
+βm,2 p (B12)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
⇒

X
A0n,2 = nm
αm,1 v21 nm
+ αm,2 v22 nm
+ βm,1 v23 nm
+ βm,2 v24 (B13)
m
 18

A01 = v11 · α1 + v12 · α2 + v13 · β1 + v14 · β2 (B14)

A02 = v21 · α1 + v22 · α2 + v23 · β1 + v24 · β2 (B15)

 
1
X V1
Ain,1 e− 2 ikn L = ×Jn−m × (B16)
m
~ω
1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (kz + χn )(kn (kz + χm ) + qm (χn − kz ) − iky (χn + χm ))
αm,1 p
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )

1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (kz + χn )(kn (kz − χm ) + qm (χn − kz ) − iky (χn − χm ))
+αm,2 p (B17)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (kz + χn )(−kn (kz + χm ) + qm (χn − kz ) + iky (χn + χm ))
+βm,1 p (B18)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (kz + χn )(kn (χm − kz ) + qm (χn − kz ) + iky (χn − χm ))
+βm,2 p (B19)
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
⇒

1
X
Ain,1 e− 2 ikn L = αm,1 unm nm nm nm
11 + αm,2 u12 + βm,1 u13 + βm,2 u14 (B20)
m

and
 
1
X V1
Ain,2 e− 2 ikn L = ×Jn−m × (B21)
m
~ω
1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (χn − kz )(−kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
αm,1 p
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )

1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (χn − kz )(kn (χm − kz ) + qm (kz + χn ) − iky (χn + χm ))
+αm,2 p (B22)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (χn − kz )(kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
+βm,1 p (B23)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (χn − kz )(kn (kz − χm ) + qm (kz + χn ) + iky (χn + χm ))
+βm,2 p (B24)
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
⇒

1
X
Ain,2 e− 2 ikn L = αm,1 unm nm nm nm
21 + αm,2 u22 + βm,1 u23 + βm,2 u24 (B25)
m

From equation (B2)

 19

 ikn L ikn L 1 1

e 2 fil11 e 2 fil21 e− 2 ikn l f0l11 e− 2 ikn L f0l21 o

Bn,1
 ik2n L ikn L 1 1  B o 
e f e 2 fil22 e− 2 ikn L f0l12 e− 2 ikn L f0l22
 ikn L il12
  n,2 
ikn L 1 1  i  =
 e 2 fil13 e 2 fil23 e− 2 ikn L f0l13 e− 2 ikn L f0l23  Bn,1
i
ikn L
e 2 fil14
ikn L
e 2 fil24
1
e− 2 ikn L f0l14
1
e− 2 ikn L f0l24 Bn,1

iqm L iqm L
 1 1

e 2 fim11 αm,1 + e 2 fim21 αm,2 + e− 2 iqm L f0m11 βm,1 + e− 2 iqm L f0m21 βm,2
iqm L iqm L 1 1
 e 2 fim12 αm,1 + e 2 fim22 αm,2 + e− 2 iqm L f0m12 βm,1 + e− 2 iqm L f0m22 βm,2
 
X  V1
 × Jn−m (B26)

 iqm L iqm L 1 1

m
 e 2 fim13 αm,1 + e 2 fim23 αm,2 + e− 2 iqm L f0m13 βm,1 + e− 2 iqm L f0m23 βm,2  ~ω
iqm L iqm L 1 1
e 2 fim14 αm,1 + e 2 fim24 αm,2 + e− 2 iqm L f0m14 βm,1 + e− 2 iqm L f0m24 βm,2

 
0
X V1
Bn,1 = Jn−m × (B27)
m
~ω
p 1
(kn + qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz + χm ) + qm (χn − kz ) − iky (χn + χm ))
αm,1 p
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )

p 1
(kn + qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz − χm ) + qm (χn − kz ) − iky (χn − χm ))
+αm,2 p (B28)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )

p 1
(kn − qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz − χn ) − iky (χn + χm ))
+βm,1 p (B29)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )

p 1
(kn − qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (χm − kz ) + qm (χn − kz ) + iky (χn − χm ))
−βm,2 p (B30)
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
⇒

X
0 nm nm nm nm
Bn,1 = αm,1 v31 + αm,2 v32 + βm,1 v33 + βm,2 v34 (B31)
m

 
0
X V1
Bn,2 = Jn−m × (B32)
m
~ω
p 1
(kn + qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (−kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
αm,1 p
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )

p 1
(kn + qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (kn (χm − kz ) + qm (kz + χn ) − iky (χn + χm ))
+αm,2 p (B33)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )

p 1
(kn − qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
−βm,1 p (B34)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )

p 1
(kn − qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz − χm ) + qm (kz + χn ) + iky (χn + χm ))
−βm,2 p (B35)
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
 20

X
0 nm nm nm nm
Bn,2 = αm,1 v41 + αm,2 v42 + βm,1 v43 + βm,2 v44 (B36)
m

B10 = v31 · α1 + v32 · α2 + v33 · β1 + v34 · β2 (B37)

B20 = v41 · α1 + v42 · α2 + v43 · β1 + v44 · β2 (B38)

 
1
X V1
i
Bn,1 e− 2 ikn L
= ×Jn−m × (B39)
m
~ω
iqm L p
(qm − kn )(qm + iky )e 2 χ3n (kz + χn )(−kn (kz + χm ) + qm (χn − kz ) − iky (χn + χm ))
αm,1 p
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )

iqm L p
(qm − kn )(qm + iky )e 2 χ3n (kz + χn )(kn (χm − kz ) + qm (χn − kz ) − iky (χn − χm ))
+αm,2 p (B40)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (kz + χn )(kn (kz + χm ) + qm (χn − kz ) + iky (χn + χm ))
+βm,1 p (B41)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (kz + χn )(kn (kz − χm ) + qm (χn − kz ) + iky (χn − χm ))
+βm,2 p (B42)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
⇒

1
X
i
Bn,1 e− 2 ikn L = αm,1 unm nm nm nm
31 + αm,2 u32 + βm,1 u33 + βm,2 u34 (B43)
m

and
 
1
X V1
i
Bn,2 e− 2 ikn L = ×Jn−m × (B44)
m
~ω
iqm L p
(qm − kn )(qm + iky )e 2 χ3n (χn − kz )(kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
αm,1 p
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )

iqm L p
(qm − kn )(qm + iky )e 2 χ3n (χn − kz )(kn (kz − χm ) + qm (kz + χn ) − iky (χn + χm ))
+αm,2 p (B45)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (χn − kz )(−kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
+βm,1 p (B46)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (χn − kz )(kn (χm − kz ) + qm (kz + χn ) + iky (χn + χm ))
+βm,2 p (B47)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
⇒
 21

1
X
i
Bn,2 e− 2 ikn L = αm,1 unm nm nm nm
41 + αm,2 u42 + βm,1 u43 + βm,2 u44 (B48)
m

From equation (B20),(B25),(B43) and (B48)

Ai1 · Mr = u11 · α1 + u12 · α2 + u13 · β1 + u14 · β2 (B49)

Ai2 · Mr = u21 · α1 + u22 · α2 + u23 · β1 + u24 · β2
B1i · Mr = u31 · α1 + u32 · α2 + u33 · β1 + u34 · β2
B2i · Mr = u41 · α1 + u42 · α2 + u43 · β1 + u44 · β2

From equation (B49)

u11 u12 u13 u14

 
u21 u22 u23 u24 
mat = 
u31 u32 u33 u34 
u41 u42 u43 u44

Inverse of the matrix becomes

γ11 γ12 γ13 γ14

 
γ21 γ22 γ23 γ24 
mat−1 =
γ31 γ32 γ33 γ34 
γ41 γ42 γ43 γ44

α1 = (γ11 · Ai1 + γ12 · Ai2 + γ13 · B1i + γ14 · B2i ) · Mr (B50)

α2 = (γ21 · Ai1 + γ22 · Ai2 + γ23 · B1i + γ24 · B2i ) · Mr
β1 = (γ31 · Ai1 + γ32 · Ai2 + γ33 · B1i + γ34 · B2i ) · Mr
β2 = (γ41 · Ai1 + γ42 · Ai2 + γ43 · B1i + γ44 · B2i ) · Mr

where

−1
γ11 = (R11 · R10 − R9−1 · R8 )−1 [R9−1 (u34 u−1 −1 −1 −1 −1
14 − R7 u24 u14 ) − R11 (u44 u14 − R12 u24 u14 )]
−1
γ12 = (R11 · R10 − R9−1 · R8 )−1 [R9−1 R7 − R11
−1
R12 ]
−1
γ13 = −(R11 · R10 − R9−1 · R8 )−1 R9−1
−1
γ14 = (R11 · R10 − R9−1 · R8 )−1 R11
−1
(B51)

−1
γ21 = (R10 · R11 − R8−1 · R9 )−1 [R8−1 (u34 u−1 −1 −1 −1 −1
14 − R7 u24 u14 ) − R10 (u44 u14 − R12 u24 u14 )]
−1
γ22 = (R10 · R11 − R8−1 · R9 )−1 [R8−1 R7 − R10
−1
R12 ]
−1
γ23 = −(R10 · R11 − R8−1 · R9 )−1 R8−1
−1
γ24 = (R10 · R11 − R8−1 · R9 )−1 R10
−1
(B52)
 22

γ31 = (R4−1 · R3 − R6−1 · R5 )−1 [R6−1 (u41 u−1 −1 −1 −1 −1

11 − R2 u21 u11 ) − R4 (u31 u11 − R1 u21 u11 )]
γ32 = (R4−1 · R3 − R6−1 · R5 )−1 [R6−1 R2 − R4−1 R1 ]
γ33 = (R4−1 · R3 − R6−1 · R5 )−1 R4−1
γ34 = −(R4−1 · R3 − R6−1 · R5 )−1 R6−1 (B53)

γ41 = (R3−1 · R4 − R5−1 · R6 )−1 [R5−1 (u41 u−1 −1 −1 −1 −1

11 − R2 u21 u11 ) − R3 (u31 u11 − R1 u21 u11 )]
γ42 = (R3−1 · R4 − R5−1 · R6 )−1 [R5−1 R2 − R3−1 R1 ]
γ43 = (R3−1 · R4 − R5−1 · R6 )−1 R3−1
γ44 = −(R3−1 · R4 − R5−1 · R6 )−1 R5−1 (B54)

R1 = (u32 − u31 u−1 −1

11 u12 )(u22 − u21 u11 u12 )
−1

R2 = (u42 − u41 u−1 −1

11 u12 )(u22 − u21 u11 u12 )
−1

R3 = u33 − u31 u−1 −1

11 u13 − R1 (u23 − u21 u11 u13 )
R4 = u34 − u31 u−1 −1
11 u14 − R1 (u24 − u21 u11 u14 )
R5 = u43 − u41 u−1 −1
11 u13 − R2 (u23 − u21 u11 u13 )
R6 = u44 − u41 u−1 −1
11 u14 − R2 (u24 − u21 u11 u14 )

R7 = (u33 − u34 u−1 −1

14 u13 )(u23 − u24 u14 u13 )
−1

R12 = (u43 − u44 u−1 −1

14 u13 )(u23 − u24 u14 u13 )
−1

R8 = u31 − u34 u−1 −1

14 u11 − R7 (u21 − u24 u14 u11 )
R9 = u32 − u34 u−1 −1
14 u12 − R7 (u22 − u24 u14 u12 )
R10 = u41 − u44 u−1 −1
14 u11 − R12 (u21 − u24 u14 u11 )
R11 = u42 − u44 u−1 −1
14 u12 − R12 (u22 − u24 u14 u12 )

recall the equation

A01 = v11 · α1 + v12 · α2 + v13 · β1 + v14 · β2 (B55)

A02 = v21 · α1 + v22 · α2 + v23 · β1 + v24 · β2
B10 = v31 · α1 + v32 · α2 + v33 · β1 + v34 · β2
B20 = v41 · α1 + v42 · α2 + v43 · β1 + v44 · β2
With the help of equation (B50), the above equation becomes

A01 =v11 · (γ11 · Ai1 + γ12 · Ai2 + γ13 · B1i + γ14 · B2i ) · Mr (B56)
+ v12 · (γ21 · Ai1 + γ22 · Ai2 + γ23 · B1i + γ24 · B2i ) · Mr
+ v13 · (γ31 · Ai1 + γ32 · Ai2 + γ33 · B1i + γ34 · B2i ) · Mr
+ v14 · (γ41 · Ai1 + γ42 · Ai2 + γ43 · B1i + γ44 · B2i ) · Mr
≡M11 · Ai1 + M12 · Ai2 + M13 · B1i + M14 · B2i
 23

Similarly

A02 = M21 · Ai1 + M22 · Ai2 + M23 · B1i + M24 · B2i (B57)
B10 = M31 · Ai1 + M32 · Ai2 + M33 · B1i + M34 · B2i
B20 = M41 · Ai1 + M42 · Ai2 + M43 · B1i + M44 · B2i
or

 i
A01 M11 M12 M13 M14 A1
  
0 i
A2  M21 M22 M23 A2i 
M24  
B 0  = M M32 M33 M34 B1 

1 31
B20 M41 M42 M43 M44 B2i
where

M11 = (v11 · γ11 + v12 · γ21 + v13 · γ31 + v14 · γ41 ) · Mr (B58)
M12 = (v11 · γ12 + v12 · γ22 + v13 · γ32 + v14 · γ42 ) · Mr (B59)
M13 = (v11 · γ13 + v12 · γ23 + v13 · γ33 + v14 · γ43 ) · Mr (B60)
M14 = (v11 · γ14 + v12 · γ24 + v13 · γ34 + v14 · γ44 ) · Mr (B61)

M21 = (v21 · γ11 + v22 · γ21 + v23 · γ31 + v24 · γ41 ) · Mr (B62)
M22 = (v21 · γ12 + v22 · γ22 + v23 · γ32 + v24 · γ42 ) · Mr (B63)
M23 = (v21 · γ13 + v22 · γ23 + v23 · γ33 + v24 · γ43 ) · Mr (B64)
M24 = (v21 · γ14 + v22 · γ24 + v23 · γ34 + v24 · γ44 ) · Mr (B65)

M31 = (v31 · γ11 + v32 · γ21 + v33 · γ31 + v34 · γ41 ) · Mr (B66)
M32 = (v31 · γ12 + v32 · γ22 + v33 · γ32 + v34 · γ42 ) · Mr (B67)
M33 = (v31 · γ13 + v32 · γ23 + v33 · γ33 + v34 · γ43 ) · Mr (B68)
M34 = (v31 · γ14 + v32 · γ24 + v33 · γ34 + v34 · γ44 ) · Mr (B69)

M41 = (v41 · γ11 + v42 · γ21 + v43 · γ31 + v44 · γ41 ) · Mr (B70)
M42 = (v41 · γ12 + v42 · γ22 + v43 · γ32 + v44 · γ42 ) · Mr (B71)
M43 = (v41 · γ13 + v42 · γ23 + v43 · γ33 + v44 · γ43 ) · Mr (B72)
M44 = (v41 · γ14 + v42 · γ24 + v43 · γ34 + v44 · γ44 ) · Mr (B73)
in general

4
X
Mij = (vil · γlj ) · Mr (B74)
l=1

For simplicity we define some matrices whose matrix elements are:

p 1
(kn − qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz − χn ) + iky (χn + χm ))
(v11 )nm = p (B75)
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )
 
V1
Jn−m
~ω
 24

p 1
(kn − qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (kz − χm ) + qm (kz − χn ) + iky (χn − χm ))
(v12 )nm = p (B76)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz + χm ) + qm (χn − kz ) + iky (χn + χm ))
(v13 )nm = p (B77)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz − χm ) + qm (χn − kz ) + iky (χn − χm ))
(v14 )nm = p (B78)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
 
V1
Jn−m
~ω

p 1
(kn − qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
(v21 )nm = − p
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )
(B79)
 
V1
Jn−m
~ω

p 1
(kn − qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz − χm ) + qm (kz + χn ) − iky (χn + χm ))
(v22 )nm = − p
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )
(B80)
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (−kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
(v23 )nm = p
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )
(B81)
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (kn (χm − kz ) + qm (kz + χn ) + iky (χn + χm ))
(v24 )nm = p (B82)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz + χm ) + qm (χn − kz ) − iky (χn + χm ))
(v31 )nm = p (B83)
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )
 
V1
Jn−m
~ω
 25

p 1
(kn + qm )(qm + iky ) χ3n (kz + χn )e− 2 iL(kn −qm ) (kn (kz − χm ) + qm (χn − kz ) − iky (χn − χm ))
(v32 )nm = p (B84)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )
 
V1
Jn−m
~ω

p 1
(kn − qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz − χn ) − iky (χn + χm ))
(v33 )nm = p (B85)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )
 
V1
Jn−m
~ω

p 1
(kn − qm )(qm − iky ) χ3n (kz + χn )e− 2 iL(kn +qm ) (kn (χm − kz ) + qm (χn − kz ) + iky (χn − χm ))
(v34 )nm = − p
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
(B86)
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (−kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
(v41 )nm = p
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )
(B87)
 
V1
Jn−m
~ω

p 1
(kn + qm )(qm + iky ) χ3n (χn − kz )e− 2 iL(kn −qm ) (kn (χm − kz ) + qm (kz + χn ) − iky (χn + χm ))
(v42 )nm = p (B88)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )
 
V1
Jn−m
~ω

p 1
(kn − qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
(v43 )nm = − p
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )
(B89)
 
V1
Jn−m
~ω

p 1
(kn − qm )(qm − iky ) χ3n (χn − kz )e− 2 iL(kn +qm ) (kn (kz − χm ) + qm (kz + χn ) + iky (χn + χm ))
(v44 )nm = − p
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
(B90)
 
V1
Jn−m
~ω

1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (kz + χn )(kn (kz + χm ) + qm (χn − kz ) − iky (χn + χm ))
(u11 )nm = p (B91)
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω
 26

1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (kz + χn )(kn (kz − χm ) + qm (χn − kz ) − iky (χn − χm ))
(u12 )nm = p (B92)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (kz + χn )(−kn (kz + χm ) + qm (χn − kz ) + iky (χn + χm ))
(u13 )nm = p (B93)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (kz + χn )(kn (χm − kz ) + qm (χn − kz ) + iky (χn − χm ))
(u14 )nm = p (B94)
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (χn − kz )(−kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
(u21 )nm = p (B95)
4kn χn (kn + iky )(qm − iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω

1 p
(kn + qm )(qm + iky )e− 2 iqm L χ3n (χn − kz )(kn (χm − kz ) + qm (kz + χn ) − iky (χn + χm ))
(u22 )nm = p (B96)
4kn χn (kn + iky )(qm − iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (χn − kz )(kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
(u23 )nm = p (B97)
4kn χn (kn + iky )(qm + iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω

iqm L p
(qm − kn )(qm − iky )e 2 χ3n (χn − kz )(kn (kz − χm ) + qm (kz + χn ) + iky (χn + χm ))
(u24 )nm = p (B98)
4kn χn (kn + iky )(qm + iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

iqm L p
(qm − kn )(qm + iky )e 2 χ3n (kz + χn )(−kn (kz + χm ) + qm (χn − kz ) − iky (χn + χm ))
(u31 )nm = p (B99)
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω
 27

iqm L p
(qm − kn )(qm + iky )e 2 χ3n (kz + χn )(kn (χm − kz ) + qm (χn − kz ) − iky (χn − χm ))
(u32 )nm = p (B100)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (kz + χn )(kn (kz + χm ) + qm (χn − kz ) + iky (χn + χm ))
(u33 )nm = p (B101)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (kz + χn )(kn (kz − χm ) + qm (χn − kz ) + iky (χn − χm ))
(u34 )nm = p (B102)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

and
iqm L p
(qm − kn )(qm + iky )e 2 χ3n (χn − kz )(kn (kz + χm ) + qm (kz + χn ) − iky (χn − χm ))
(u41 )nm = p (B103)
4kn χn (kn − iky )(qm − iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω

iqm L p
(qm − kn )(qm + iky )e 2 χ3n (χn − kz )(kn (kz − χm ) + qm (kz + χn ) − iky (χn + χm ))
(u42 )nm = p (B104)
4kn χn (kn − iky )(qm − iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (χn − kz )(−kn (kz + χm ) + qm (kz + χn ) + iky (χn − χm ))
(u43 )nm = p (B105)
4kn χn (kn − iky )(qm + iky ) χ3m (kz + χm )
 
V1
× Jn−m
~ω

1 p
(kn + qm )(qm − iky )e− 2 iqm L χ3n (χn − kz )(kn (χm − kz ) + qm (kz + χn ) + iky (χn + χm ))
(u44 )nm = p (B106)
4kn χn (kn − iky )(qm + iky ) χ3m (χm − kz )
 
V1
× Jn−m
~ω

1
(Mr )nm = e− 2 ikn L δn, m (B107)

[1] U. Fano, Phys. Rev. 124, 1866 (1961).

 28

[2] C.-C. Lee and C.-L. Ho, Annals of Physics 320, 175 (2005).
[3] S.-L. Zhu and Z. D. Wang, Phys. Rev. B 65, 155313 (2002).
[4] J.-H. Dai and R. Zhu, The European Physical Journal B 87, 288 (2014).
[5] R. Zhu, J.-H. Dai, and Y. Guo, Journal of Applied Physics 117, 164306 (2015).
[6] R. Zhu and M. Lai, Journal of Physics: Condensed Matter 23, 455302 (2011).
[7] R. Zhu and J. Berakdar, Phys. Rev. B 81, 014403 (2010).
[8] M. Moskalets and M. Büttiker, Phys. Rev. B 70, 245305 (2004).
[9] W. Li and L. E. Reichl, Phys. Rev. B 60, 15732 (1999).
[10] E. Tekman and P. F. Bagwell, Phys. Rev. B 48, 2553 (1993).
[11] J. Cui, B. Ji, X. Song, and J. Lin, Plasmonics 14, 41 (2019).
[12] R. Zhu and C. Cai, Journal of Applied Physics 122, 124302 (2017).
[13] B. Sun, L. Zhao, C. Wang, X. Yi, Z. Liu, G. Wang, and J. Li, The Journal of Physical Chemistry C 118, 25124 (2014).
[14] S. Karmakar, R. K. Varshney, and D. R. Chowdhury, OSA Continuum 2, 531 (2019).