Diffusion Monte Carlo calculations of the polarizability of a confined hydrogen atom: benchmarking and application to high symmetry wells

We present a non-perturbative direct method to calculate the polarizability of a hydrogen atom confined in a three-dimensional potential well of any geometry. The calculation is based on the diffusion Monte Carlo method. The advantage of the method is simplicity of implementation and immediate adaptability to any well shape. The method is validated for the well-studied spherically confined hydrogen atom, and demonstrated in the case of two other geometries (cube and octahedron), for which this paper provides the first set of results. Although demonstrated here for the confined hydrogen atom, the method can be immediately applied to any single-electron system placed in a three-dimensional potential well of any type and geometry. Results for a hydrogen atom confined in potential wells of cubic or spherical symmetry suggests that the polarizability in these cases is a universal function of the volume of the well. This result can simplify calculations for real situations such as in quantum dots.


Introduction
Polarizability is a fundamental quantity characterizing atomic and molecular systems, and one of the first examples of quantum calculation performed on such systems (Landau and Lifshitz 1977).It is connected to optical properties and polarization interactions (Chandrasekhar 1960).Many past studies have shown that atomic and molecular polarizability is strongly reduced by confining the system into a potential well (Kirkwood 1932, Buckingham 1937, Fowler 1984, Holka et al 2005).This effect, relevant from a fundamental point of view as well as for its applications, has been studied with different methods applied to several systems, in particular the confined hydrogen atom (Banerjee et al 2002, Montgomery andSen 2012).This system was originally introduced to simulate the effect of high pressure on its static dipole polarizability (Michels et al 1937).
Most of the studies about the dipole polarizability of confined hydrogen atom make use of the approximate Kirkwood formula (Kirkwood 1932), where the dipole polarizability is estimated from the radial expectation value r 2 á ñ of the unperturbed wave functions of the atom in consideration.Actually, it has been demonstrated (Dalgarno 1962) that, for most atomic systems in S states, the Kirkwood polarizability underestimates the exact value .
a Calculations of the polarizability a of a compressed hydrogen atom have been reported in (Dutt et al 2001), using the approximate variational wave functions and the Kirkwood (Dalgarno 1962) and Buckingham (Buckingham 1937) approximations. (Sen et al 2002) report calculations of static dipole polarizability for a shellconfined hydrogen atom, using uncoupled perturbation-numerical calculations.Static and dynamic dipole polarizabilities, calculated using variational perturbation theory, can be found in (Montgomery 2002).
In this paper, we present a method for the direct calculation of the polarizability of a confined hydrogen atom: a finite but weak external electric field is applied to the confined system, and the resulting dipole moment is calculated using the diffusion Monte Carlo (DMC) method.We report the results of the direct calculations of the polarizability a of a hydrogen atom confined in boxes of different shapes.Our technique is easily implemented and, if requested, is easily extended to calculate the components of the polarizability tensor for asymmetric systems and nonlinear terms of the polarizability.It is validated here for a spherical well and demonstrated for a cubic and an octahedral well.

Method
The framework of the method is the determination of the polarizability from the induced dipole moment of the quantum system by a small, but finite applied electric field.This procedure is not new, but our proposal makes it very usable and versatile.The dipole moment is calculated from the wave function, that comes from a DMC calculation and is sampled on a three-dimensional mesh.The mesh used is the simplest Cartesian one, regardless of the actual geometry of the system under consideration, and this allows immediate extensibility to the nontrivial tensor case and to the non-linear case.The resulting method has also the advantage of not requiring any set of basic functions.
The method employs the possibilities offered by the DMC algorithm to calculate exactly the wave function of a one-or two-electron (in a singlet state) system confined in a potential well of any shape, based on the formal replacement of the Schrödinger equation with a reaction-diffusion equation in an imaginary time coordinate.
The direct calculation implies knowledge of the charge density as a function of the three spatial coordinates.The dipole moment is given by the integral: where the charge density is in atomic units, and x y z , , ( ) y comes from the sampling of the walker's positions on the mesh, for a long enough imaginary time span, as described below.
For systems whose polarizability tensor is diagonal, such as those that will be considered in this work, the dipole moment d must be parallel to the electric field .
E This was verified in all cases presented here (see next section).The scalar polarizability α is then given by (Sukumar and Navaneethakrishnan 1989): As it is known, DMC realizes the reaction-diffusion process of numerous (typically, 10 3 ~) 'walkers' in the phase space: a walker is a mathematical point in a set {R i } of phase coordinates R in a 3N-dimensional space, where N is the number of electrons in the system.Walkers diffuse stochastically based on the short-time diffusion-reaction propagator (Foulkes et al 2001): where t is a numerical step.The first exponential term represents the quantum kinetic energy terms, while the second exponential term defines the birth/death algorithm applied to the walkers.The estimated energy E T controls the total population of the walkers, and it is adapted by a relaxation method until the numerical equilibrium of the population of walkers is reached.In this way, the positions of the walkers at steady state sample the wave function r .

( ) y
To increase simplicity and adaptability, DMC is used here in its unbiased version, described in the authors' previous works (Micca Longo et al 2015, Micca Longo et al 2021a, Micca Longo et al 2021b), which avoids the common use of a priori functions (importance sampling).Since such functions need to be reformulated for any new geometry, the unbiased version allows to simply change the well geometry and, if requested, the number and position of atomic nuclei in the wells.To speed up calculations, a different formulation of the short time propagator (4) is used, equivalent in the diffusion regime: this is obtained replacing the diffusion factor in equation (4) by a uniform distribution in a cube centered into r and with edge length 2 3 , t to match the variance of the gaussian propagator (Micca Longo et al 2021a).
As previously mentioned, since in DMC the walkers sample the wave function, not the charge density, a numerical grid in three-dimensional space is needed to sample the number of walker steps during a sufficiently long simulation, after the steady state of the reaction and diffusion kinetics is reached.The square of the function is then calculated from the sampled values.The same grid is used to numerically calculate the integral expression of d.The electric field is specified by its magnitude E and by the set of direction cosines.
The contribution of the applied electric field E to the Hamiltonian H operator used in the DMC calculations is given by r.
• E - The full Hamiltonian is therefore where V is the potential energy term, and it is equal to or , +¥ depending on r being inside or outside the well.Note that the first term of H is not actually requested by the DMC method, since its effect on y is obtained by the formal diffusion of walkers.
The DMC calculation is non-perturbative and indeed exact: this means that our method allows us to calculate all the non-linear components of the polarizability.To determine polarizability in the strict sense, a small enough electric field ( 0.02 E = atomic units, a. u.) is applied so that nonlinear effects are negligible.The field magnitude cannot be reduced arbitrarily, though, since the signal/noise ration of the stochastic determination of the polarization is negatively affected.The value of the polarizability is determined within its second significant figure: this precision, however, should be more than sufficient for a value of polarizability for use in concrete calculations, for example of optical properties.

Benchmarking
In this section, we show that the proposed method reproduces the results obtained in the literature in the case of spherical confinement.In this case, V in equation ( 5) can be written explicitly as where c is the confinement radius.
Figure 1 refers to a spherical confinement and shows the polarization calculated as a function of the radius of the spherical well.Both quantities are in atomic units.This calculation uses 5 10 3 ´walkers, for a total imaginary time of 10 6 a. u., with a calculation step of 10 2 -a.u. for larger confinements and of 10 3 -for tighter ones.
The intensity of the applied electric field is a u 0.02 .., E = and we have verified that the results are not modified by a substantial reduction of this value, even if this reduction involves an equally substantial increase in the statistical fluctuations of the resulting dipole moment.
The mesh for the sampling of x y z , , ( ) y covers a cubic domain embedding the well, using 40 ´40 ´40 elements.This number has been optimized through numerical tests: increasing the number of mesh elements does not have a significant effect on the sampled quantities, while it increases statistical fluctuations thus requiring more calculation time.Results have been validated with a comparison with (Cohen et al 2008, Saha et al 2011, He et al 2023), finding a quite good agreement, as shown in figure 1.The residual difference in the case of very large confinement radii is relatively unimportant, since we are interested in the behavior produced by the confinement, that is reproduced very accurately.
A residual difference is observed only for large confinement radii, for which the atom is essentially free.This residual difference appears here since the mesh used does not accurately resolve the peak of the wave function near the origin, in a free atom.On the other hand, a much finer mesh leads to aforementioned disadvantage.One solution could be the use of an adaptive mesh.Although this can be verified in future work, it does not appear to be necessary for the present discussion, since the behavior produced by the confinement is accurately reproduced.
This result should also be evaluated in the light of the simplicity of implementation of the proposed method, in contrast to the heavily analytical methods used in the benchmark calculations.In this regard, we emphasize that the numerical code implementing our approach only needs to realize the reaction-diffusion process described in the previous section.It is a very streamlined code, relatively easy to write and requiring no special functions or coordinates.Of course, this approach is not aimed at specifying numerous decimal places of the result, which is often not necessary in concrete applications.

Application to high symmetry wells: hints of universal behavior
In this section, new results are given for two cases for which no data is available in the literature: the cubic confinement and the octahedral confinement.Both cases are potentially important in determining physical properties of hydrogen atoms trapped in crystal lattices (Pupyshev 2011, de Graaf et al 2020, Gend et al 2023).Both cube and octahedron have a O h symmetry, and this leads to scalar polarizability.In the case of cubic and octahedral confinement, the horizontal axis is defined by the explicit expression of V in the two cases, respectively: Correspondingly, in the case of the cube, 2c is the edge length; in the case of the octahedron; 2c is the edge of the circumscribed cube.Since the choice of these lengths for the two solids is somewhat arbitrary, it is convenient to use the volume of the polyhedron as a parameter, which is always well defined.This parameter is convenient also because polarizability is dimensionally a volume.As we will see later, this choice brings additional consequences.Figure 2 shows the polarizability of a hydrogen atom confined in a cube and in an octahedron, respectively.There is an explicit demonstration concerning O h symmetry of both cases: for any electric field direction, the dipole moment stays parallel to the field and its magnitude does not change.Here, the electric field is directed along symmetry axes (C 4 , C 3 , C 2 ).The possibility of changing the electric field direction as desired, with respect to the potential well, is obviously a useful feature in other less symmetrical cases.Indeed, the proposed direct method can be applied to very different confinement geometries, as demonstrated in our previous works where well shapes as different as cylindrical and icosahedral were considered (Micca Longo et al 2021a, Micca Longo et al 2021b, Longo et al 2023).Furthermore, in the case of a confinement for which the system has a non-trivial polarizability tensor, the calculation of its different components is straightforward.
The method presented in this paper allows us to highlight an issue that may be of significant value in applying the concept of an atom confined in a potential well to real systems.As mentioned above, the linear dimension of a potential well of generic shape can be defined in several ways (such as diameter, radius, edge length, diagonal length).This requires some attention when examining a graph of the polarizability of a quantum system as a function of a linear dimension.Accordingly, in figure 2, we use well volume as a parameter.For a sufficiently symmetric well, such as T d of O h symmetry, the volume completely specifies the well; this is not the case, for example, of a cylindrical geometry.
This choice makes the three cases comparable and we performed this comparison.The resulting graph is shown in figure 3: bearing in mind the noise inherent in this type of simulations, the three lines using the same abscissa defined as above are almost coincident.At the very least, they are quite close, within the first decimal place in atomic units, sufficiently precise for many applications of polarizability.
The results of our calculations thus show that the polarizability of a hydrogen atom confined in a potential hole of sufficiently symmetrical shape may depend only on the volume of the well, or at least is strictly parameterized by it.This observation, that probably deserve to be substantiated by future strictly theoretical or analytical considerations, may be useful in applying the concept to real systems.For example, in the study of an atomic impurity at a quantum dot of sufficiently symmetric geometry, the polarizability of the system can be calculated in the spherical case and then related to the same volume of the system considered.

Conclusions
In this paper, an approach to directly calculate the polarizability, or the components of the polarizability tensor, is proposed, in the case of a one-electron system, and specifically a hydrogen atom, confined in a threedimensional potential well.The approach is simple to implement and adapt, and it has the advantage of using Cartesian coordinates in any case, not requiring any expansion in basic functions.The atomic nucleus does not  necessarily have to be placed in the center of the well.Results show that the polarizability of a hydrogen atom confined in an infinitely deep and sufficiently symmetric potential well may depend primarily on its size and not on its shape.This is a potentially useful result for applications to nanotechnology, since it allows to use results of a spherical well to estimate those of a highly symmetric but non spherical well of the same volume.This result was made possible thanks to the flexibility of our approach to any shape of the well, with minimal changes to the code that implements it.Although this has not been demonstrated in this work, the approach can also be applied to confined hydrogen-like atomic ions and even to a molecular system such as H 2 + .For these different applications, the only modification required is to rewrite the potential part V of the Hamiltonian, affecting a single line of code.

Figure 1 .
Figure 1.Polarizability a of a hydrogen atom confined in a sphere of radius c: comparison between literature (Cohen et al 2008, Saha et al 2011, He et al 2023) data sets and DMC results.

Figure 2 .
Figure 2. Polarizability , a as a function of the volume V, of a hydrogen atom confined in a cubical (left) and in an octahedral well (right).The label of the figure refers to O h symmetry axes (C 4 , C 3 and C 2 ) chosen as the direction of the electric field.

Figure 3 .
Figure 3. Calculated polarizability a of the hydrogen atom confined into potential wells of three different geometries, plotted as a function of the well volume.This plot evidences the substantial coincidence of the three results (see text for discussion).