# Semi-classical limit for the varying-mass Schrödinger equation with random inhomogeneities

The varying-mass Schrödinger equation (VMSE) has been successfully applied to model electronic properties of semiconductor hetero-stuctures, for example, quantum dots and quantum wells. In this paper, we consider VMSE with small random heterogeneities, and derive radiative transfer equations for its solutions. The main tool is to systematically apply the Wigner transform in the semiclassical regime (the rescaled Planck constant ε≪ 1), and then expand the resulted Wigner equation to proper orders of ε. As a proof of concept, we numerically compute both VMSE and radiative transfer equations, and show that their solutions agree well.

## Authors

• 14 publications
• 38 publications
• 29 publications
12/26/2020

### Semi-classical limit of an inverse problem for the Schrödinger equation

It is a classical derivation that the Wigner equation, derived from the ...
03/29/2021

### Higher dimensional generalization of the Benjamin-Ono equation: 2D case

We consider a higher-dimensional version of the Benjamin-Ono (HBO) equat...
10/07/2019

### A particle method for the homogeneous Landau equation

We propose a novel deterministic particle method to numerically approxim...
07/10/2021

### Unconditionally stable exponential time differencing schemes for the mass-conserving Allen-Cahn equation with nonlocal and local effects

It is well known that the classic Allen-Cahn equation satisfies the maxi...
12/03/2019

### Transforming the Lindblad Equation into a System of Linear Equations: Performance Optimization and Parallelization

Rapidly growing performance and memory capacity of modern supercomputers...
11/28/2021

### Remarks on the Radiative Transfer Equations for Climatology

Using theoretical and numerical arguments we discuss some of the commonl...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

The Schrödinger equation with varying mass has gained great attention in solid state physics, and been successfully used to to study electronic properties of semiconductor hetero-structures [27, 24, 11, 29]

. For example, it can describe localized defects in crystalline media, which may yield bound states localized to the defect. It is also related to describe non-compact line defect (edge) perturbations, e.g., by interpolating two-dimensional honeycomb structures via domain walls

[13, 12]. Edge modes may be produced by such perturbations, which propagate in the direction parallel to the edge and localize transverse to the edge.

We are interested in deriving the asymptotic limit of the following varying-mass Schrödinger equation

 iε∂tuε(t,x)+12ε2∇x⋅(m(t,x)∇xuε(t,x))=0, (1.1)

where , with and is the rescale Planck constant. The varying mass , which can be time-dependent, is assume to be random and highly oscillatory, with a given covariance matrix in time and space. We shall assume that decays fast enough at infinity to validate all the derivations. One goal of the paper is to show that in the regime, the Wigner transform of the solution converges to a special radiative transfer equation.

The problem is motivated by a fact that simulating (1.1) is extremely challenging in the semi-classical regime (). The challenges are two-folded. In deterministic regime (meaning is a deterministic highly oscillatory function in ), standard numerical solvers require to the small wavelength of to be resolved, for example, a mesh size and time step of order is required when finite difference methods are used [22, 23]. The time-splitting spectral method [4, 6] can improve the mesh size to be of order , however, has limitations to compute the Schrödinger equation with varying mass. A bigger problem comes from the randomness in . Since only the covariance of

is given, numerically one has to find many realizations and compute the deterministic Schrödinger equation before finding the ensemble mean/variance of the solution. The number of realizations, however, increases as

, as details in the random fluctuation become more and more important.

Disregarding the challenges from the stochasticity, merely for the deterministic system, alternative approaches have been developed. These included the WKB-type methods, e.g., Gaussian beam methods [18, 17] and frozen Gaussian approximation [15, 19, 20]. The idea is to apply the WKB-type ansatz

 uε(t,x)=A(t,x)exp(iS(t,x)ε),

and derive the eikonal equation for and transport-like equation for , where both and are functions of large scale, i.e., -independent. To our best knowledge, no such types of methods have been applied to efficiently solve (2.1) in the literature yet. And even for standard Schrödinger equation with random potential term, the application of the methods have not been fully understood.

In the paper, we shall systematically derive asymptotic equations for (2.1) by the Wigner transform [14], which is a main tool in semiclassical theory parallel to the WKB-type methods mentioned above. The literature on deriving the asymptotic equations for wave propagation in random media [2, 3, 1, 26, 21, 10, 9, 8, 7] is rich, most of which started with the Schrödinger equation with constant mass, and the randomness and high oscillations are introduced through the potential term. When it is the effective mass term that is random and highly oscillatory, the process of the derivation is rather similar but much more delicate, as will be detailed later in our paper. As a proof of concept, we numerically verify the derived radiative transfer equation by carefully computing and comparing its solution to the one of VMSE (1.1), and show that the two solutions agree.

The rest of the paper is organized as follows. To better illustrate the derivation, we start with a simpler case with being deterministic and only spatially dependent, and derive the limiting radiative transfer equation by the Wigner transform in Section 2. In Section 3, we systematically introduce the derivation of the limiting equation for the varying-mass Schrödinger equation (1.1) with random heterogeneities. We present our numerical validation in Section 4 and make conclusive remarks in Section 5.

## 2 Wigner transform of VMSE in the deterministic setting

As a preparation, we first investigate the semi-classical limit for (1.1) with deterministic mass in this section, which is to consider:

 iε∂tuε(t,x)+12ε2∇x⋅(m0(x)∇xuε(t,x))=0. (2.1)

The varying mass is a real function of . It is assumed to be deterministic and time-independent. We study the Cauchy problem, i.e., the is assumed to decay at infinite. is a complex function and we typically care about its physical observables such as the energy density and the energy flux , defined respectively as

 ρε(t,x)=|uε(t,x)|2,Jε(t,x)=εIm(¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯uε(t,x)∇xuε(t,x)).

Wigner transform is a technique explored in [25] for the Schrödinger equation with random potential, and has been demonstrated as a very powerful tool for investigating the semi-classical limit. It defines a function on the phase space:

 Wε(t,x,k)=1(2π)d∫Rdeikyuε(t,x−ε2y)¯¯¯¯¯uε(t,x+ε2y)dy. (2.2)

Here is the complex conjugate of

. This definition is essentially the Fourier transform of

 ⟨x−ε2y|u⟩⟨u|x+ε2y⟩

on the variable.

The Wigner transform loses phase information: if is perturbed to

, the Wigner transform is kept the same. However, the physical observables can be recovered, namely, the first and second moments of

provide the energy density and the energy flux:

 ∫RdWε(t,x,k)dk=ρε(t,x),∫RdkWε(t,x,k)dk=Jε(t,x). (2.3)

It is not guaranteed that is positive, and thus it does not serve directly as the energy density on the phase space. By plugging in the Schrödinger equation, we derives the equation satisfied by in the following Lemma.

Let satisfy the VMSE (2.1), then its Wigner transform (2.2) satisfies:

 ∂tWε+1εQε1Wε+Qε2Wε=εQε3Wε, (2.4)

where the three operators are defined to be:

 Qε1Wε=|k|22∫Rdeipx(2π)d~m1(t,p)i[Wε(t,x,k−ε2p)−Wε(t,x,k+ε2p)]dp (2.5)
 Qε2Wε=k2⋅∫Rdeipx(2π)d~m1(t,p)[∇xWε(t,x,k−ε2p)+∇xWε(t,x,k+ε2p)]dp (2.6)
 Qε3Wε= 18∫Rdeipx(2π)d~m1(t,p)i[ΔxWε(t,x,k−ε2p)−ΔxWε(t,x,k+ε2p)]dp (2.7) +18∫Rdeipx(2π)d~m1(t,p)i|p|2[Wε(t,x,k−ε2p)−Wε(t,x,k+ε2p)]dp.

where the Fourier transform of in space is defined by

 ~m0(t,p)=∫Rde−ipzm0(t,z)dz. (2.8)

The proof is direct derivation. Notice that

 ∂tWε=1(2π)d∫Rdeiky∂tuε¯¯¯¯¯uεdy+1(2π)d∫Rdeikyuε∂t¯¯¯¯¯uεdy, (2.9)

we have, plugging in (2.1):

 ∂tWε= iε2(2π)d∫Rdeiky∇x⋅(m0(x−ε2y)∇xuε(x−ε2y))¯¯¯¯¯uε(x+ε2y)dy (2.10) −iε2(2π)d∫Rdeiky∇x⋅(m0(x+ε2y)∇x¯¯¯¯¯uε(x+ε2y))uε(x−ε2y)dy := iε2(2π)dM1−iε2(2π)dM2.

Since the two terms and are conjugate with for the second term, we only study the first one. With integration by parts:

 M1= 2ε∫Rd[∇y(eiky)⋅∇xuε(x−ε2y)]m0(x−ε2y)¯¯¯¯¯uε(x+ε2y)dy (2.11) +2ε∫Rdeikym0(x−ε2y)∇xuε(x−ε2y)⋅∇y¯¯¯¯¯uε(x+ε2y)dy := I1+I2.

We treat the and respectively in the following. Perform integration by parts again to

 I1= 4ε2∫RdΔy(eiky)m0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy (2.12) +4ε2∫Rd∇y(eiky)⋅∇ym0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy +2ε∫Rd[∇y(eiky)⋅∇x¯¯¯¯¯uε(x+ε2y)]m0(x−ε2y)uε(x−ε2y)dy := I11+I12+I13.

Note that and the last term can be combined so that a complete gradient of is available, namely one arrives at a formula for

 I1= 12I1+12(I11+I12+I13)=12(I1+I13)+12(I11+I12) (2.13) = 2ε2∫RdΔy(eiky)m0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy + 2ε2∫Rd∇yeiky⋅∇ym0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy + 1ε∫Rd∇y(eiky)⋅∇x[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]m0(x−ε2y)dy.

For in (2.11), integration by parts against produces

 I2= 4ε2∫Rd[∇y(eiky)⋅∇y¯¯¯¯¯uε(x+ε2y)]m0(x−ε2y)uε(x−ε2y)dy (2.14) + 4ε2∫Rdeiky[∇ym0(x−ε2y)⋅∇y¯¯¯¯¯uε(x+ε2y)]uε(x−ε2y)dy + ∫Rdeikym0(x−ε2y)uε(x−ε2y)Δx¯¯¯¯¯uε(x+ε2y)dy:=I21+I22+I23.

On the other hand, integration by parts against gives

 I2= 4ε2∫Rd[∇y(eiky)⋅∇yuε(x−ε2y)]m0(x−ε2y)¯¯¯¯¯uε(x+ε2y)dy (2.15) + 4ε2∫Rdeiky[∇ym0(x−ε2y)⋅∇yuε(x−ε2y)]¯¯¯¯¯uε(x+ε2y)dy + ∫Rdeikym0(x−ε2y)¯¯¯¯¯uε(x+ε2y)Δxuε(x−ε2y)dy:=I′21+I′22+I′23.

Note that and can be combined after another integration by parts

 I21+I′21= −4ε2∫Rd∇y(eiky)⋅∇ym0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy (2.16) − 4ε2∫Rdm0(x−ε2y)Δy(eiky)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy.

and can be combined similarly

 I22+I′22= −4ε2∫Rd∇y(eiky)⋅∇ym0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy (2.17) −

Hence using (2.14)-(2.17) and the trick in (2.13), one derives the formula for in (2.11)

 I2= 14(I21+I22+I23)+12I2+14(I′21+I′22+I′23) (2.18) = (14I23+12I2+14I′23)+14(I21+I′21)+14(I22+I′22) = 14∫Rdeikym0(x−ε2y)Δx[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy − 2ε2∫Rd∇y(eiky)⋅∇ym0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy − 1ε2∫Rdm0(x−ε2y)Δy(eiky)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy −

Finally from (2.13) and (2.18), one gets

 M1= 14∫Rdeikym0(x−ε2y)Δx[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy (2.19) +1ε∫Rdm0(x−ε2y)∇y(eiky)⋅∇x[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy +1ε2∫Rdm0(x−ε2y)Δy(eiky)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy −1ε2∫RdeikyΔym0(x−ε2y)[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy := T1+T2+T3+T4.

All the terms can be explicitly expressed by the Wigner transform (2.2). In particular:

 T1 =∫Rdeipx~m0(p)ΔxWε(x,k−ε2p)dp, (2.20) T2 =∫Rdeipx~m0(p)ik⋅∇xWε(x,k−ε2p)dp, T3 =∫Rd−|k|2eipx~m0(p)Wε(x,k−ε2p)dp, T4 =ε2∫Rd|p|2eipx~m0(p)Wε(x,k−ε2p)dp.

We use as an example to show this. Recalling:

 ΔxWε(x,k)=1(2π)d∫eikyΔx[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy,

we have

 ∫Rdeipx~m0(p)ΔxWε(x,k−ε2p)dp (2.21) = 1(2π)d∫∫∫eipxe−ipzm0(z)ei(k−ε2p)yΔx[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dzdpdy, = ∫Rdeikym0(x−ε2y)Δx[uε(x−ε2y)¯¯¯¯¯uε(x+ε2y)]dy=14T1,

where we used the fact that

 δ(x)=1(2π)d∫Rdeixzdz,and1(2π)d∫∫f(x)eixzdzdx=f(0). (2.22)

Using (2.20), we get

 M1= 14∫Rdeipx~m0(p)ΔxWε(x,k−ε2p)dp (2.23) +1ε∫Rdeipx~m0(p)ik⋅∇xWε(x,k−ε2p)dp +∫Rd|p|2eipx~m0(p)Wε(x,k−ε2p)dp

By the conjugate argument, one gets, setting :

 M2= 14∫Rdeipx~m0(p)ΔxWε(x,k+ε2p)dp (2.24) −1ε∫Rdeipx~m0(p)ik⋅∇xWε(x,k+ε2p)dp +∫Rd|p|2eipx~m0(p)Wε(x,k+ε2p)dp

Finally, substitute (2.23) and (2.24) into (2.10), and we arrive at the Wigner equation in (2.4).

[Formal] With certain regularity, the formal limit of (2.4) is a transport equation:

 ∂tW+m0(x)k⋅∇xW−|k|22∇xm0(x)⋅∇kW+O(ε2)=0. (2.25)

The trajectory of particles follows:

 ˙x=km0(t,x),˙k=−|k|22∇xm0(t,x), (2.26)

and

 |k|2m0(x)=const. (2.27)

To prove the limiting equation (2.25), one only needs to derive the limiting behavior of . Indeed, formally, assuming :

 1εQε1Wε=−|k|22∫Rdeipx(2π)d~m1(t,p)ip⋅∇kWε(t,x,k)dp+O(ε2)=−|k|22∇xm0⋅∇kWε+O(ε2),

and that

 Qε2Wε=k⋅∫Rdeipx(2π)d~m1(t,p)∇xWε(t,x,k)dp+O(ε2)=m0k⋅∇xWε+O(ε2).

One then arrives at (2.25) by plugging them in the original equation (2.4). To show (2.27) with , one simply takes the time derivative along the trajectory:

 ddt[|k(t)|2m0(x(t))]=2k⋅˙km0+|k|2∇xm0⋅˙x=0,

according to the trajectory equation (2.26), and thus the quantity is conserved along the trajectory.

This proposition essentially guarantees the positivity of the solution to the limiting equation (2.25).

## 3 Semi-classical limit for VMSE with random perturbation

We consider the VMSE where the effective mass involves random perturbation, namely:

 iε∂tuε+12ε2∇x⋅(mε(t,x)∇xuε)=0, (3.1)

where the effective mass is

 mε(t,x)=m0(t,x)+√εm1(t/ε,x/ε). (3.2)

While the leading order is assumed to be deterministic and smooth, we allow the random perturbation to present small scales at . Furthermore we assume it is mean-zero and stationary in both and with the correlation function :

 R(t,x)=E[m1(s,z)m1(t+s,x+z)]∀x,z∈Rd and t,s∈R. (3.3)

Taking the Fourier transform of the function in both time and space, one has:

 ^R(ω,p)=∫Rd+1e−iωs−ipzR(s,z)dsdz, (3.4)

then it is straightforward to show:

 E[~m1(τ,p)^m1(ω,q)]=(2π)de−iωτ^R(ω,p)δ(p+q), (3.5)

and

 ^R(−ω,p)=^R(ω,p),and^R(ω,−p)=^R(ω,p).

We dedicate this section to the derivation of the semi-classical limit of the equation above. We will show that In the zero limit of , the Wigner transform of , the solution to the VMSE (2.1) with varying random mass (3.2), solves the radiative transfer equation:

 ∂tW+m0k⋅∇xW−k22∇xm0⋅∇kW=1(2π)d∫Rd14(p⋅k)2^R(m02(p2−k2),p−k)[W(p)−W(k)]dp. (3.6)

In view of (2.4) in Lemma 2, noting that , the Wigner equation (3.1) is transformed to

 ∂tWε+1εQε1Wε+Qε2Wε+1√εPε1Wε+√εPε2Wε=εQε3Wε+1√εPε3Wε, (3.7)

where the operators are defined in (2.5)-(2.7), and are their counterparts defined by :

 Pε1Wε=|k|22∫Rdeipξ(2π)d~m0(τ,p)i[Wε(k−12p)−Wε(k+12p)]dp, (3.8)
 Pε2Wε=k2⋅∫Rdeipξ(2π)d~m0(τ,p)[∇xWε(k−12p)+∇xWε(k+12p)]dp, (3.9)

and

 Pε3Wε= ε28∫Rdeipξ(2π)d~m0(τ,p)i[