Time-Domain Simulation of Three Dimensional Quantum Wires

A method is presented to calculate the eigenenergies and eigenfunctions of quantum wires. This is a true three-dimensional method based on a direct implementation of the time-dependent Schrödinger equation. It makes no approximations to the Schrödinger equation other than the finite-difference approximation of the space and time derivatives. The accuracy of our method is tested by comparing it to analytical results in a cylindrical wire.


Introduction
Quantum wires have become extremely important structures in the field of quantum nanodevices. They have found applications in such diverse fields as quantum electronics [1][2][3][4][5] and nonlinear optics [6][7][8][9][10]. This paper describes the use of the finite-difference time-domain (FDTD) method to determine the eigenenergies and eigenfunctions of quantum wires from the potential function [11][12][13][14][15]. Specifically, a method is described to accurately model long, thin wires without using excessive computer resources. Furthermore, a Bessel function method is presented that can be used to verify the calculations of the FDTD method.
An important application of this method is to nonlinear optics, where nonlinear susceptibilities are calculated from the quantum properties of the material. For example, transition moments, which are calculated from wave functions, and energies are used in a sum-overstates expression to get nonlinear susceptibilities. This approach was used to find optimized potentials [16,17] as well as to analyze the character of a potential that gives the largest nonlinearity by analyzing ones made of a finite number for straight line segments, as demonstrated by Burke [18][19][20].
FDTD can be used to quickly calculate nonlinear susceptibilities, which can be used to optimize the nonlinear response by varying the shapes of arbitrary structures. In this work, we consider structures made of fused conducting cylinders because these are the types of structures that can be made using common fabrication techniques found in nanotechnology labs. As such, this work can be seen as providing the computational building blocks with applications to nonlinear optics, electronics, and related fields.
In the next section of this paper, we describe the implementation of the FDTD method and its use to calculate eigenenergies and eigenfunctions. This includes a description of a method to increase the resolution of the FDTD simulation by averaging the potential in cells bordering two different materials. We then present the results of determining the eigenenergies of a quantum wire at various radii along with the Bessel function results to confirm the accuracy. The flexibility of the method is demonstrated by modeling potentials that do not otherwise easily lend themselves to analysis. We end by summarizing the results of the paper and discussing future directions. S1 Appendix describes the Bessel function method.

The Finite-Difference Time-Domain Method
We begin describing the FDTD implementation of the time-dependent Schrödinger equation [13,14], which is written in the following form @cðx; y; z; tÞ @t ¼ i ℏ 2m e @ 2 cðx; y; z; tÞ @x 2 þ @ 2 cðx; y; z; tÞ @y 2 þ @ 2 cðx; y; z; tÞ @z 2 À i ℏ Vðx; y; zÞcðx; y; z; tÞ: We separate ψ(x,y,z,t) into real and imaginary components: cðx; y; z; tÞ ¼ c real ðx; y; z; tÞ þ i Á c imag ðx; y; z; tÞ: Inserting Eq 2 into Eq 1 leads to two coupled equations: @c real ðx; y; z; tÞ @t ¼ À ℏ 2m e r 2 c imag ðx; y; z; tÞ þ 1 ℏ Vðx; y; zÞc imag ðx; y; z; tÞ; ð3aÞ @c imag ðx; y; z; tÞ @t ¼ ℏ 2m e r 2 c real ðx; y; z; tÞ À 1 ℏ Vðx; y; zÞc real ðx; y; z; tÞ: To code these equations, we take the finite-difference approximations in space and time which results in the following two coupled equations: c kþ1 real ðm; n; lÞ ¼ c k real ðm; n; lÞ þ Dt ℏ Vðm; n; lÞc kþ1=2 imag ðm; n; lÞ À ℏ 2m e Dt ðDxÞ 2 c kþ1=2 imag ðm þ 1; n; lÞ þ c kþ1=2 imag ðm À 1; n; lÞ þc kþ1=2 imag ðm; n þ 1; lÞ þ c kþ1=2 imag ðm; n À 1; lÞ þc kþ1=2 imag ðm; n; l þ 1Þ þ c kþ1=2 imag ðm; n; l À 1Þ À 6c kþ1=2 imag ðm; n; lÞ ; ð4aÞ c kþ3=2 imag ðm; n; lÞ ¼ c kþ1=2 imag ðm; n; lÞ À Dt ℏ Vðm; n; lÞc kþ1 real ðm; n; lÞ Dt ðDxÞ 2 c kþ1 real ðm þ 1; n; lÞ þ c kþ1 real ðm À 1; n; lÞ þc kþ1 real ðm; n þ 1; lÞ þ c kþ1 real ðm; n À 1; lÞ þc kþ1 real ðm; n; l þ 1Þ þ c kþ1 real ðm; n; l þ 1Þ À 6c kþ1 real ðm; n; lÞ In Eq 4 integer indices m,n.l representing the positions in a matrix have replaced the Cartesian coordinates x,y,z, respectively, in Eq 3. Similarly, the time step k has replaced t. Once the cell size Δx is chosen, the time step Δt must also be chosen so the constants preceding the spatial Laplacian are small enough to maintain stability [21]. The alternate iteration of the real and imaginary Eq (4a) and (4b) simulates the behavior of the waveform propagating in time. Details are available in the literature [11][12][13][14][15]. It is also necessary to be able to absorb outgoing waveforms to prevent them from being reflected back into the problem space and interfering with the simulation. This is accomplished with a perfectly matched layer (PML) [22]. The PML effectively adds artificial impedance between the real and imaginary parts of the waveform. The PML being used is based on Zheng [23] as described in a previous paper [24], so details will not be repeated here. Fig 1 illustrates the quantum wires to be simulated. The aspect ratio is at least ten to one, so we will use a length of 100 Angstroms, and a diameter of 10 Angstroms or less. We assume that the inside of the wire is at zero potential and the background medium is 4.6 eV corresponding to the work function of silver. Fig 2 illustrates the problem space being used to simulate quantum wires. Normally, the potential of a quantum structure is specified by setting the V parameter in Eq 4a and 4b) to the potential corresponding to the location in the structure being simulated. As an example, the cross-section of a quantum wire of radius 5 Angstroms using cells of one Angstrom might appear as shown in Fig 3A. The red indicates the inside of the wire which is at zero potential, while the green indicates the surrounding potential at V = 4.6 eV. Grey represents the PML. Borrowing from a method developed for electromagnetic simulation [25], a more accurate representation of the boundaries can be achieved by averaging the fraction of each material that is in a cell and assigning the value of potential accordingly. The results are shown in Fig 3B.

Determining the Eigenenergies and Eigenfunctions
The FDTD method can be used to determine the eigenenergies and eigenstates of a potential that does not otherwise lend itself to an analytic solution. Any quantum wave function in a given quantum system can be written in the following manner, where the ϕ n (r) are the eigenfunctions of the system and the ε n are the corresponding eigenenergies [11]. We can write the state variable in this fashion even if we do not know the eigenfunctions and eigenenergies. The eigenenergies can be determined by monitoring the time-domain data at one point in the problems space, say r 0 , and then taking the Fourier transform, The last step results from the following Equations: The above transform produces a series of delta functions in the frequency domain corresponding to the eigenenergies.
The simulation starts by initializing a test function in the potential, as shown in Fig 4. The test function is a narrow Gaussian pulse. It chosen to ensure that it contains all the eigenfunctions of the structure being analyses. As the time-domain simulation proceeds, the pulse spreads out, some of it being absorbed by the PML. However, most of the waveform remains in the wire. During the simulation, the time-domain data at the original source point is saved. All the waveform simulations in this paper are in three dimensions. Since only two dimensions are displayed, we show the middle of the problems space in the y direction.
After 100,000 iterations, the program is halted. The stored time-domain data is shown in Fig 5(A), and the Fourier transform of this data is shown in Fig 5(B). The first two eigenenergies appear at 0.6378 and 0.6488. In order to verify the accuracy of the FDTD technique, an analytic method based on Bessel functions has been developed and is described in S1 Appendix of this paper. Table 1 shows a comparison of the FDTD determination of the ground state energies of the quantum wire of To determine the eigenfunction ϕ m (r) corresponding to an eigenenergy ε m , we take the discrete Fourier transform of the state variable at the frequency ω m = ε m / ℏ at every point in the problems space:  To construct the ground state eigenfunction at ε 0 = 0.6378 eV the original test function is initialized in the problem space as in

Other Wire Configurations
In the previous section the simulation examples used simple cylinders to establish the accuracy by comparison with an analytic solution. The true strength of the FDTD method is that it can be used on a potential that is not necessarily a simple, regular structure. Two examples are given in this section.
By adding short cylinders perpendicular to the main cylinder, it may be possible to substantially increase the hyperpolarization of the system [6]. One possible configuration is illustrated in Fig 8(A), where a small cylinder seven Angstroms long was added to the main cylinder. Fig 8  (B) shows the corresponding problem space. Fig 8C illustrates the corresponding ground state eigenfunction.
Another example is shown in Fig 9(A) where the cylinder is tapered slightly in the middle. The tapering has been proposed as a method of controlling transport in nanowire MOSFETs (metal oxide semiconductor, field effect transistors.) [2].

Conclusion
We have shown that the FDTD method can be used to determine the eigenenergies and eigenfunctions of quantum wires on size scales that are of interest to nanotechnologists. We have applied the FDTD method to solve the simple cylinder, which has analytic solutions in the form of Bessel functions. The two methods give eigenenergies that are accurate to better than two percent. Comparisons with other methods were not made. The more complex structures that we have calculated here cannot be tested against analytic solutions. However, given that they are made of merged cylinders, we expect that the accuracy will be comparable.
The simulations described in this paper were done on an HP DL140 with eight cores, a high end workstation, and does not represent extraordinary computational resources. A typical simulation to determine the eigenenergy and eigenfunction requires about 40 seconds. Note that these simulations require no assumptions about the shape of the object, so computational speeds are not affected by the shape.
Future projects include simulating multiple particles in a quantum wire. The FDTD method has already been shown to be capable of simulating two particles in a quantum dot [26].