1.1. Potential wells with simple geometry

Often a real-world object can be modeled with some exceptionally simple physical system. An electron in a potential well with the shape of an exact cube or a sphere, and the potential energy V stepping abruptly at the well boundary is one example. In that situation, the potential energy of the electron is (spherical case)

or (cubic case)

Here the V ₀ is the depth of the well. Θ is the Heaviside theta function,

ρ is the radius of the circular or spherical potential well and L is the side length of the cubic potential well. The variable r is just the distance from the origin in 2D or 3D space: r=x2+y2 or r=x2+y2+z2. The function V appears in the time-independent Schrödinger equation, from which the energy eigenvalues and eigenfunctions are calculated

Here x is a two- or three-dimensional vector.

1.2. Lamé circles and spheres

In real-world situations it can not always be assumed that the electron is confined to an exactly spherical or cubic space. The choice V ₀ = ∞ is another crude approximation. A logical improvement to the model would be to define 2D curves and 3D surfaces with a shape between a square and a circle or cube and a sphere. One geometrical object of this type is the Lamé circle (or supercircle)

and the Lamé sphere (or supersphere)

Yet another name for these curves and surfaces is superquadrics. The objects described by these equations, with α = b = c = 1 (as in all cases studied in this article), become a circle or sphere when q = 2 and a square or a cube when q → ∞. The parameter r > 0 is the radius of the circle or sphere that these shapes become when q → 2. As a general result, an n-dimensional solid of this shape has volume [⁹]

where Γ is the gamma function. To produce a set of surfaces described by (6), having different values of q but same constant enclosed volume V, the parameter r in Eq. (6) has to be made a function of q:

or the equivalent for a 2D Lamé circle.

3D surfaces of shapes given by Eq. (6) with α = b = c = 1 are shown in Fig. 1.

FIGURE 1 Wolfram Mathematica 11.2.0 plots of the surface defined by |x| ^q + |y| ^q + |z| ^q = r ^q for values q = 0.5, q = 1.0 q = 1.5, q = 2.5, q = 3 and q = 4 from left to right and up to down. The parameter r has been set to value 1. 

The system with a potential well bounded by a Lamé curve has been studied earlier by Bera et al. [¹⁰], but only for the 2D case and with infinite V ₀. In this article, the results will be extended to the equivalent 3D systems with several finite values of V ₀ allowed.

A similar method of analytical geometry exists for converting a sphere continuously into a Platonic solid (regular octahedron, dodecahedron or icosahedron instead of a cube) [¹¹,¹²], and it can be expected that the energy eigenvalues of an electron in that kind of potential well approach the spherical case when more faces are added to the polyhedron.

1.3. Experimental methods

A practical system that can be modelled as one or more electrons in a 1 to 3-dimensional potential well is the quantum dot (QD), which is a nanometer-scale semiconductor particle having an absorption/emission spectrum similar to a system of electrons confined in a potential well of equivalent size. An electron in a large QD has a lower ground state energy than one in a small QD, just like an electron in a short versus long 1D square well. One complication in this theory of QDs is the -possibly position-dependent- difference of electron and hole effective masses to their normal value or the value in bulk semiconductor [¹³]. This kind of nanoparticles can be produced by chemical or mechanical means such as molecular beam epitaxy or chemical vapor deposition [¹⁴], and their shape can be probed by methods such as measuring their spin and parity oscillations in a magnetic field [¹⁵]. Information about the effect of shape changes, as in the model systems of this article, on the electronic spectrum of the system can also be utilized for quantum dot shape determination. Conversely, this also means that QDs with desired spectral properties can be produced by tuning their shape in an appropriate way.

1.4. Computational approach

The methods applied in this article are numerical, with ground state and first excited state energies obtained with quantum diffusion Monte Carlo (QDMC) or Implicit finite differencing (IFD) on imaginary time axis. The system is a 2D or 3D potential well with boundary described by Eqs. (5) and (6). The parameters α, b and c all have value 1 , so the well is not oblong in any particular direction. The system contains one or two confined electrons. A suitable nonlinear fitting function is sought for the data points of ground state energy E ₀ as a function of q and V ₀. This is done to make it easier to estimate the E ₀ for other values of q, V ₀ by interpolation. With spectroscopic data and these results, it may be possible to estimate whether an electron-trapping cavity in a real material is more close to sphere or cube in shape.

Both the QDMC and IFD methods for finding the ground state are based on the time evolution of an initial state |Ψ(0)⟩= c ₀|ψ ₀⟩ + c ₁|ψ ₁⟩ + c ₂|ψ ₂⟩... written as a linear combination of the eigenstates of H^, |ψ _n ⟩:

Here it is assumed that H^ has no explicit time dependence in Schrödinger picture. Defining an imaginary time coordinate with equation t = −iτ, the time evolution in terms of τ is

Provided that the E _n form an increasing sequence of positive numbers, the state vector |Ψ(τ)⟩ of (10) approaches something proportional to |ψ ₀⟩ when τ → ∞. This allows the numerical calculation of the ground state wave function ψ ₀(x) = ⟨x| ψ ₀⟩. The eigenvalue E ₀ can be obtained by following the expectation value ⟨Ψ(τ)|H^|Ψ(τ)⟩/⟨Ψ(τ)|Ψ(τ)⟩, or the value of H^Ψ(x,τ)/Ψ(x,τ) at some constant position x, during the simulation.

The difference between the QDMC and IFD methods is that in QDMC the diffusive imaginary-time dynamics is modelled as a system of random walkers with a source term allowing creation and destruction of them. The IFD method calculates actual wave functions that do not have to be real and positive valued.

2.1. Transcendental equation and shooting method for q = 2 and q = ∞

2.2 Quantum DMC for one particle in 2D and 3D

The diffusion Monte Carlo results for both 2D and 3D were calculated with values of 2 ≤ q ≤ 20 and 50 ≤ V ₀ ≤ 20000.

The result for each combination of q and V ₀ was calculated 4 times and averaged to one data point. The time step length in both 2D and 3D simulations was ∆t = 10⁻⁵, the energy update parameter was α = 40 and the maximum allowed number of random walkers was 10⁵ (2D) or 9 × 10⁵ (3D). Simulated time length was 20 units, and 1/4 from the beginning of the energy output file was discarded before averaging the remaining energy values to an estimate for E ₀.

The results were plotted as graphs where either V ₀ or q was kept constant and the other varied. Particle mass was set to 1 in all calculations. The area of a 2D potential well was kept constant at A = 4 when changing the value of the exponent q, and the volume of the 3D potential well was kept at V = 8. This was done by defining r = r(q) as in (8).

Calculations for well depths of 600, 2000 and 20000 were also done similarly for a 3D well of volume V=2.0 to see the effect of volume change on the ground state energy. In the V=2.0 calculation, the time step was ∆t = 5 × 10⁻⁶, the maximum number of walkers was 90000, length of simulation was 9.0 time units and the energy update parameter α had value 100.0.

In preliminary calculations, it was seen that an exponential function of type

fits well in the data points of E ₀ as a function of q, and a rational square root function

to the data points of E ₀ as a function of V ₀. Making Eq. (20) more complex by converting it to a sum of two exponentials did not add any more accuracy in the fit. Both kind of graphs (as a function of q or V ₀) can be used for estimating the ground state energy for intermediate values of q and V ₀ by interpolation. There is no obvious theoretical reason for the E ₀(V ₀) and E ₀(q) to be close to these functions; the reason for calculating the nonlinear fit is just to predict the behavior of the functions between calculated data points.

The results were plotted and the iterative nonlinear fits calculated with the Linux application Grace-5.1.25. The obtained values of fitting parameters A ₀, A ₁ and A ₂ were recorded for all cases. Another way to do the curve fit is to fix A ₀ and A ₁ to values predicted by a transcendental equation or with the shooting method, and leave only A ₂ free.

2.3. Quantum DMC for two particles in 2D

2.4. Implicit finite-differencing for one particle in 2D

3.1. Ground state energy of q = 2 and q = ∞ potential wells

The results for the E ₀ of square/cubic systems with side length 2, obtained from the graphical solution of Eq. (12), are on the left side of Table I.

The equivalent results for 2D circular or 3D spherical potential wells, with area 4 or volume 8 are shown on the right side of Table I. The former has been obtained with the shooting method and the latter with Eq. (16).

3.2. QDMC results for one particle in 2D or 3D

For both 2D and 3D, and for every value of V ₀ from 50 to 20000, 12 ground state energy data points from q = 2 to q = 20 were obtained in the QDMC calculation. To each point set of different V ₀, the exponential function of (20) was fitted to the E ₀(q) data. This was done both with all parameters A ₀, A ₁ and A ₂ free, or with A ₀ and A ₁ fixed to the limiting values of E ₀(q = 2) and E ₀(q → ∞) for that V ₀ (Table I) and only A ₂ left free. Figures 2 and 3 show both the E ₀ data points and the 1- or 3-parameter curve fits to the data with a potential well of area A = 4.0 or volume V=8.0.

FIGURE 2 QDMC ground state energies and 1- or 3-parameter exponential curve fits (20) for one electron in a 2D potential well of shape (5) and area 4.0, and with different values of exponent q. The units have been chosen so that ℏ=me=1 and energy is dimensionless. 

FIGURE 3 QDMC ground state energies and 1- or 3-parameter exponential curve fits (20) for one electron in a 3D potential well of shape (6) and volume 8.0, and with different values of exponent q. The dimensionless units are again ℏ=me=1. 

The similarity of the 1-parameter (gray lines) and 3parameter (black lines) exponential fits in Figs. 2 and 3 reflects how accurately the E ₀ values calculated with QDMC agree with those calculated by different means for the q = 2 and q = ∞ limiting cases (Table I).

Main important features in the data of Fig. 2 and 3 include that the E ₀ grows monotonically with increasing q in all cases, and the difference in E ₀ between q = 2 and q = ∞ is in the range of 5-15%. Most of the change in ground state energy in all data sets has already taken place when q has increased to value 8 from the original 2. The nonlinear fit parameter A ₂ was between −0.30 and −0.20 for almost all data sets.

Figure 4 shows the additional results E ₀(q) for the three cases of a 3D potential well of smaller volume 2.0. The graphs look similar to those for larger potential wells, except for the E ₀ values being higher. The curve fit parameters of Figs. 2 to 6 are also included as tables of numerical data in the Appendix.

FIGURE 4 QDMC E ₀ values and 3-parameter exponential curve fits (20) for one electron in a 3D potential well of shape (6) and volume 2.0, and with exponent q as the independent variable. 

It is also possible to present the data as graphs with q and V kept constant and V ₀ varied. The square root type function (21) was fitted iteratively to E ₀(V ₀) data sets corresponding to different values of q. The parameters A ₀ and A ₁ in Eq. (21) were free, and not pre-fixed, in each case. Figures 5 and 6 show the data points and curve fits, with V ₀ on the horizontal axis and E ₀ on the vertical.

FIGURE 5 QDMC ground state energies of one electron as function of V ₀ for a supercircular potential well of area 4.0, shape given by Eq. (5) and several values of exponent q. The square root type curve fits of Eq. (21) are also shown in the figure. 

FIGURE 6 QDMC ground state energies of one electron as function of V ₀ for a superspherical potential well of volume 8.0, shape of Eq. (6) and several values of exponent q. The square root type curve fits of Eq. (21) are also shown in the figure. 

3.3. QDMC results for two particles in 2D

The QDMC ground state energies, E ₀, for two interacting or non-interacting electrons in a 2D potential well of depth V ₀ = 200 or V ₀ = 20000 are shown in Fig. 7 with the value of parameter s varying from 1.0 to 4.0. The area of the potential well is A = π when s = 1 regardless of what the value of q is and A = πs ² for other values of s, to make the results easier to compare. For q = 2 it holds that s = r, for q = 8 the relation is s ≈ 1.1162r and for q = ∞, s ≈ 1.1284r with r as in Eq. 5. More generally, for a given value of q and 2p-dimensional space the connection between s and r is s=Vol(r,2,q)/π with Vol(r,n,q) as in Eq. (7).

FIGURE 7 QDMC ground state energy E ₀ for two interacting or non-interacting (NI) electrons in a supercircular potential well of depth V ₀ = 200 (left) or V ₀ = 20000 (right). The parameter q is the same as in Eq. (5) and s a linear scaling factor, with the area of the supercircle being πs ² for a given value of s. An E ₀(s) ∝ s ⁻² nonlinear fit has been made to the q = 8 data points to demonstrate the difference caused by the presence of Coulomb interaction. 

In the graphs, there is also a curve fit of form E ₀(r) ∝ r ⁻² made for each q = 8 data point set. The energy eigenvalues of a set of non-interacting particles, in a potential well deep enough for approximation V ₀ = ∞ to be valid, are proportional to reciprocal square of the well diameter. From Fig. 7 it is easy to see that this type of function can be quite accurately fitted to the q = 8 data points for noninteracting electrons, even for the lower well depth V ₀ = 200, but not accurately at all to the point sets with Coulomb potential included.

The expectation values hr ₁₂i are shown as function of s for potential wells of both depths in Fig. 7. Parabolas ⟨r ₁₂⟩ (s) = c ₂ s ² + c ₁ s + c ₀ have also been fitted in each set of data points, and the result is that each point set is close to being on the same line with |c ₁| > 50|c ₂|.

The effect of making the parameter ( (meant to prevent division by zero when calculating the electron-electron Coulomb potential) smaller by a factor of 100 was also tested in the q = 2, V ₀ = 200 case with two interacting electrons. Figure 8 shows how much difference this causes in the E ₀(s) data points. The value ( = 0.001 may look large, especially when it’s inside the square root in r12=(x1-x2)2+(y1-y2)2+ϵ, but the graph in Fig. 8 shows that not much improvement can be made by setting a smaller value.

FIGURE 8 QDMC expectation values hr ₁₂i for two interacting or non-interacting (NI) electrons in a supercircular potential well of depth V ₀ = 200 (left) or V ₀ = 20000 (right), parameter q as in Eq. (5) and s a linear scaling factor such that the area of the supercircle is πs ² for a given value of s. 

The Appendix contains the tables of numerical data for E ₀ (Fig. 7) and ⟨r ₁₂⟩ (Fig. 8) as a function of r and for different well depths.

3.4. Implicit time stepping results for one particle in 2D

When exponential functions of type (20) were fitted in the E ₀ and E ₁ values obtained with the IFD method and imaginary time propagation, the parameters A ₀, A ₁ and A ₂ of Eq. (20) were found to be as in Table II. The data points E ₀(q) and E ₁(q) are shown with the nonlinear fits in Fig. 10.

FIGURE 9 The values of E ₀(s), obtained with QDMC, for two interacting electrons in a potential well with q = 2 and V ₀ = 200 and with two different values of cutoff parameter (. 

FIGURE 10 Energy eigenvalues E ₀ (left) and E ₁ (right) for one electron in a 2D potential well of shape (5), area 4.0 and depth V ₀ = 2000 or V ₀ = 20000. The results were calculated by imaginary time propagation done with implicit finite differencing and operator splitting. In the left image, the equivalent QDMC results are also plotted for comparison. 

The A ₀ values (limit of E at q → ∞) for the ground states of V ₀ = 2000 and V ₀ = 20000 are 2.3896 and 2.4366. When these are compared to values in Table I, 2.3912 and 2.4430, the error in the first case is less than 0.1% and that in the second, less than 0.3%. The corresponding ground state A ₁ values are 2.2072 and 2.2448, representing the limit of E ₀ for q → 2. The equivalent values in Table I, calculated with the shooting method, are 2.2087 and 2.2511. The difference between the IFD and shooting method values is less than 0.3% in both.

The A ₀ values for the first excited state, calculated with IFD method, were 5.9779 and 6.0957. The 1st excited state energy E ₁ can also be calculated for the q → ∞ system with Eqs. (12) and (13), and the results for V ₀ = 2000 and V ₀ = 20000 are 5.9779 and 6.1073. The differences between the IFD results and those from (12) and (13) are here less than 0.2%.

The value for the second excited state energy E ₂ from the IFD simulations done for q = 2 case were 11.6338 for V ₀ = 2000 and 11.8318 for V ₀ = 20000. When a shooting method calculation was done for these two systems with a rotation invariant 2nd excited state, the results were 11.635 for V ₀ = 2000 and 11.861 for V ₀ = 20000. Here the percentual error is again less than 0.3%. The consistency between q = 2 and q = ∞ results obtained with different methods makes it seem likely that the results of Fig. 10 are close to correct also for values of q between the two extremes.