Optimal Gaussian basis set for hydrogen atom in magnetic field

Question

Copy from here https://mathematica.stackexchange.com/questions/285044/optimal-basis-set-of-gaussian-functions-for-describing-a-quantum-system

This question arose during the discussion of the previous questions: https://mathematica.stackexchange.com/questions/284809/why-does-the-minimum-eigenvalue-change-dramatically-when-one-basis-function-is-a#284841 and https://mathematica.stackexchange.com/questions/284910/finding-excited-states-using-the-condition-of-wave-functions-orthogonality

Brief description of the problem:

Consider the following system with the Hamiltonian: $H=-\frac{1}{2}\Delta-\frac{1}{r}+\frac{25}{8}\rho^2-5/2$, where $r=(\rho,z,\phi)$ is a coordinate in the cylindrical system. Physically, this Hamiltonian describes a hydrogen atom in a magnetic field equal to 5 (in dimensionless units).
In order to solve this problem by the matrix method, it is necessary to choose a set of basis functions. To describe ground and excited states with the z component of the angular momentum m = 0 (1s, 2s, 3d (m = 0), 3s ...), a Gaussian basis set is well suited. Such basis functions have the following form: $\psi_j=e^{-b_{j} z^2}e^{-a_{j} \rho^2}$, where $a_{j}$ and $b_{j}$ are parameters. The task is reduced to the correct finding of parameters.

The parameters describing the ground state (1s) with high accuracy were found in this work (see here https://arxiv.org/abs/1709.05553 or here https://pubs.aip.org/aip/jcp/article/147/24/244108/195534/Accurate-and-balanced-anisotropic-Gaussian-type). But the basis presented in this article does not describe excited states (2s, 3d (m = 0), 3s ...) well, since there are not enough basis functions to describe them. The basis set parameters presented in this article has a complex structure, so increasing the basis set functions in this case is a difficult task.

It seems that a more general way to find the basis set of parameters is presented in this article (see here https://drive.google.com/file/d/1q74hAn0UAdNd8DtPkoCsdYPkr61-xhr2/view?usp=sharing or here https://onlinelibrary.wiley.com/doi/abs/10.1002/pssb.2220930140). The sets of parameters that the authors use are geometric progressions:
$\alpha_j=\alpha_{j-1}(\frac{\alpha_N}{\alpha_1})^{1/(N-1)}$

$\beta_j=\alpha_{j-1}(\frac{\beta_M}{\beta_1})^{1/(M-1)}$,
where $\alpha_1$, $\alpha_N$, $\beta_1$ and $\beta_M$ are the first and the last parameters. About finding these parameters, the authors write the following: "the parameters $\alpha_1$, $\alpha_N$, $\beta_1$ and $\beta_M$ were varied in a four dimensional optimization search which minimized the sum of the eigen-values of interest" and "We have performed the calculations presented in this paper with a basis with ten Gaussians in the $\rho$ direction and twelve Gaiissians in the z direction, for a total of 120 different basis functions.". I really don't understand what that means?

Energy values from the articles:
ground state (1s) = - 1.3793 (exact value = -1.380398866427)
the first excited state (2s) = -0.19335 (exact value = -0.193746709717)
the second excited state (3d (m=0)) = -0.07365
the third excited state (3s) = -0.03835

All these states should have the same structure of the basis set function ($\psi_j=e^{-b_{j} z^2}e^{-a_{j} \rho^2}$) due to the m=0 (the z component of the angular momentum).

I have a few questions:
I really can't understand what the authors mean in varied optimization search of the parameters $\alpha_1$, $\alpha_N$, $\beta_1$ and $\beta_M$, could you please explain to me what that means?

Why do the parameters ($\alpha_1$, $\alpha_N$, $\beta_1$ and $\beta_M$) have different indexes, shouldn't their number be the same? Since they are paired in the basis function.

I would be glad to any comments on solving the problem.

I would like to draw attention to the fact that to get the energy from the Aldrich, C., & Greene, R. L article it's needed the value from the table 1 to subtract the field value and divide by 2 (this is due to the fact that the authors of this article do not take into account the spin). For example for 2s state: value from the table 1 for magnetic field = 5 is equal 4.6133 so the energy of the 2s state = (4.6133 - 5)/2 = -0,19335

Answer 1

I really can't understand what the authors mean in varied optimization search of the parameters α1, αN, β1 and βM, could you please explain to me what that means?

In doi:10.1002/pssb.2220930140 Aldrich and Greene write "For given quantum numbers $m$ and $q$, the parameters $\alpha_1$, $\alpha_N$, $\beta_1$, and $\beta_N$ were varied in a four dimensional optimization search which minimized the sum of the eigenvalues of interest."

This means that they minimize the sum of the eigenvalues by changing the values of the four parameters appropriately. They do not mention what algorithm they use for the optimization (it's an old article), but you could do this nowadays with e.g. gradient descent methods or the Nelder-Mead "amoeba" method.

Why do the parameters (α1, αN, β1 and βM) have different indexes, shouldn't their number be the same? Since they are paired in the basis function.

They are four different parameters. The ansatz they use for the exponents is known in quantum chemistry as even-tempering [1]: the $i$th exponent is given by $$\alpha_{i} = \alpha_{0} \beta^{i}$$ when one indexes $i\in[0,N-1]$ as in C or Python, or $$\alpha_{i} = \alpha_{1} \beta^{i-1}$$ when using $i \in [1,N]$ as in Fortran. This trivially yields the relation $$\beta = \left( \frac {\alpha_{N}} {\alpha_{1}} \right)^{1/(N-1)} $$ that you have copied above.

The only difference to this case is that the magnetic field acts like a confining potential in directions orthogonal to the magnetic field [2], which means that they employ two sequences of exponents: one for the functions along $z$, the direction of the magnetic field axis, and another along $\rho$ which is the polar radial coordinate in $(x,y)$ that are the directions orthogonal to the field.

"We have performed the calculations presented in this paper with a basis with ten Gaussians in the ρ direction and twelve Gaiissians in the z direction, for a total of 120 different basis functions.". I really don't understand what that means?

Aldrich and Greene use a tensorial basis set in $\rho$ and $z$, see their equation (4): $$ \chi_\alpha (\rho,z) = \chi^{\rho}(\rho) \chi^{z}(z).$$ In addition, you also have the $\phi$ direction; however, since the Hamiltonian does not depend on this variable, Noether's theorem leads to it being a constant of motion and the $\phi$ solution is of the form $\exp(im\phi)$ as also happens for linear molecules. [3]

The $\rho$ direction is represented with 10 basis functions, $\{\chi^{\rho}_\mu\}_{\mu=1}^{10}$. The $z$ direction is represented with 12 basis functions, $\{\chi^{z}_\mu\}_{\mu=1}^{12}$. Forming all products, you get 120 basis functions.

I would like to draw attention to the fact that to get the energy from the Aldrich, C., & Greene, R. L article it's needed the value from the table 1 to subtract the field value and divide by 2 (this is due to the fact that the authors of this article do not take into account the spin)

Their calculations do take into account the spin; the spin operator even appears in the Hamiltonian [2]! The division by two is just a consequence of unit conversion: they report energies in Rydberg, and $ 1 \text{Ry} = 0.5 E_h$.

References

1. See for instance Theor. Chim. Acta 52, 231 (1979)
2. Lehtola et al, Mol. Phys. 118, e159798 (2020)
3. S. Lehtola, Int. J. Quantum Chem. 119, e25968 (2019)