TutorialGST01_GS_SE

Tutorial 1: Schrödinger-equation ground state

This tutorial walks through the process of setting up and calculating the Schrödinger-equation ground state of a one-dimensional hydrogen atom model.

Introduction

The Schrödinger-equation ground-state corresponds to the lowest energy solution to the eigen-value problem

$$\hat{\mathcal{H}} \psi (x)=E\psi (x),~~~~~~(1)$$

where $\hat{\mathcal{H}}$ is the Schrödinger-equation Hamiltonian operator, $\psi$ is the wave function, and $E$ its associated energy. In atomic units, the Hamiltonian operator is

$$\hat{\mathcal{H}} =-\frac{\Delta }{2}+\hat{\mathcal{V}} .~~~~~~(2)$$

where $-\frac{\Delta }{2}$ is the kinetic operator and $\hat{\mathcal{V}}$ is the atomic or molecular potential operator.

In this tutorial, we illustrate how to use the QMol-grid package to calculate the ground-state wave function of a one-dimensional hydrogen-like atom. Specifically, it walks through defining (i) the domain and grid discretization over which the Schrödinger equation and wave function are calculated, the (ii) atomic potential and (iii) Schrödinger-equation model, and calculating (iv) the ground state associated with these properties.

Setting up the Schrödinger-equation system

Here we present a minimal example for modeling an atomic system. We discuss additional options and more advanced features in the next tutorial.

Potential

We model the 1D hydrogen model atom using a soft-Coulomb potential $V(x)=-1/\sqrt{x^2 +a^2 }$ [Javanainen 1988] with

H = QMol_Va_softCoulomb('softeningParameter',sqrt(2));

where we choose the softening parameter $a=\sqrt{2}$ to match H's ground state energy. By default, the atom is located at the origin $x=0$ .

Note that H only corresponds to the atomic model, which is shared with molecular systems and various quantum frameworks. Thus, it must be turned into a valid Schrödinger-equation potential, using

V = QMol_SE_V('atom',H);

Domain

The simulation domain must be a Cartesian grid -- with all increasing, equally spaced discretization points -- and should be wide enough and with small enough of a discretization step to properly capture the wave function. In our case, we select a domain ranging -15 to 15 a.u., with a discretization step of 0.1 a.u.

x = -15:.1:15;

Schrödinger-equation model

We now have all the elements to define a Schrödinger-equation model object with the potential and domain defined above

SE = QMol_SE(                        ...
         'xspan',                x,  ...
         'potential',            V);

To check the configuration of our Schrödinger-equation object, we use the run-time documentation feature; this requires for the SE object to be initialized first

SE.initialize;
SE.showDocumentation;

yielding (removing the license and funding disclaimers)

=== External components ==================================================
  * convertUnit
    V-01.02.000 (09/30/2022)                                     F. Mauger
  * fourierTool
    V-01.02.002 (07/18/2023)                                     F. Mauger

=== Theory ===============================================================
  * Schrodinger equation (SE)
    V-01.21.000 (06/17/2024)                                     F. Mauger
  * The variational SE model matches the canonical Hamiltonian formalism
    describing the time evolution of the system [Mauger 2024].
    V-01.21.000 (06/17/2024)                                     F. Mauger

=== Discretization =======================================================
  * Domain discretization                                   Cartesian grid
    Grid = -15:0.1:15
    Size = 301 (7 x 43) points
    V-01.21.000 (06/17/2024)                                     F. Mauger

=== Potential ============================================================
  * Potential
    V-01.21.000 (06/17/2024)                                     F. Mauger
  * Atom-like center(s)
    > ???, parameterized as                                 (soft Coulomb)
      Z =  1.00 | a =  1.41 | X0 =   0.00
  * Soft-Coulomb potential [Javanainen 1988]                (soft Coulomb)
    Parameterized as V(x) = -Z ./ sqrt( (x-X0).^2 + a^2 ). 
    V-01.21.000 (06/17/2024)                                     F. Mauger

=== System model =========================================================
  * Electronic structure                                    wave functions
    1 wave function(s)

=== References ===========================================================
  [Javanainen 1988] J. Javanainen, J.H. Eberly, and Q. Su, "Numerical
    simulations of multiphoton ionization and above-threshold electron
    spectra," Physical Review A 38, 3430 (1988).
  [Mauger 2024] F. Mauger, C. Chandre, M.B. Gaarde, K. Lopata, and K.J. 
    Schafer, "Hamiltonian  formulation and symplectic split-operator 
    schemes for time-dependent  density-functional-theory equations of 
    electron dynamics in molecules," Communications in Nonlinear Science 
    and Numerical Simulation 129, 107685 (2024).
  [Mauger 2024b] F. Mauger and C. Chandre, "QMol-grid: A MATLAB package
    for quantum-mechanical simulations in atomic and molecular systems," 
    arXiv:2406.17938 (2024).

The first two blocks ("External components" and "Theory") list the components used in the SE objects. For most users, it can be used to (i) check the level of theory and (ii) get relevant references associated with the specific flavor of simulations used in the model.

The third block ("Discretization") summarizes the domain discretization, including the prime decomposition for the number of grid points. This decomposition can be used to optimize the fast-Fourier transforms used to evaluate the Laplacian and other spatial-derivative operators.

The fourth block ("Potential") summarizes the model potential. Molecular models contains a list of atomic centers, which are listed first and then followed by a description of the different types of pseudopotentials. In our case we have a single soft-Coulomb atomic center. The "???" results from not giving a name to the atomic model H and can be ignored.

The fifth block ("System model") summarizes the electronic structure for the Schrödinger-equation model. A Schrödinger-equation object can support several wave functions, all sharing the same potential, like various electronic states in an atom or molecule. The default number of wave functions is one.

The last block ("References") contains the list of references cited in the run-time documentation.

Ground-State Calculation

With the Schrödinger-equation object defined above, we next move to calculating its associated ground-state wave function and energy using the two commands

GSS = QMol_SE_eigs;
GSS.computeGroundState(SE);

The first command creates a Schrödinger-equation eigen-state solver GSS, which is then used on our system and prints out (again, removing the disclaimers)

=== External components ==================================================
  * convertUnit
    V-01.02.000 (09/30/2022)                                     F. Mauger
  * fourierTool
    V-01.02.002 (07/18/2023)                                     F. Mauger

=== Build Schrodinger-equation (SE) model ================================
  * Discretization                                                      OK
  * Wave function(s)                                                    OK
  * Potential                                                           OK

  * Eigen-state solver for SE Hamiltonians            MATLAB eigs function
    Tolerance  = 1e-12 
    Max. iter. = 300
    Basis dim. = 100
    V-01.21.001 (07/01/2024)                                     F. Mauger

=== References ===========================================================
  [Mauger 2024b] F. Mauger and C. Chandre, "QMol-grid: A MATLAB package
    for quantum-mechanical simulations in atomic and molecular systems," 
    arXiv:2406.17938 (2024).

=== Funding ==============================================================
    The original development of the QMol-grid toolbox, and its (TD)DFT
  features, was supported by the U.S. Department of Energy, Office of
  Science, Office of Basic Energy Sciences, under Award No. DE-SC0012462.
    Addition of the (TD)SE features was supported by the National Science
  Foundation under Grant No. PHY-2207656.

=== Wave-function energies ===============================================
  Wave fcn      Energy (-eV)         Error(a.u.)
  --------     ------------          -----------
      1           13.606              3.213e-13
  ----------------------------------------------

=== Schrodinger-equation-component energies ==============================
  Component      Energy (a.u.)      Energy (eV)
  -----------    -------------     -------------
  Kinetic             0.073              1.999
  Potential          -0.573            -15.605
  -----------    -------------     -------------
  Total              -0.500            -13.606
  ----------------------------------------------

Like in the run-time documentation, the first block ("External components") specifies the versions for some of the components used in the ground-state calculation.

All ground-state calculations start by (re)initializing the Schrödinger-equation object. The second block shows the status for this (re)initialization as well as the parameters used by the ground-state solver.

The last-but-one block ("Wave-function energies") displays the energy of the eigen state, as well as an estimate of the eigen-state error. This estimate is defined as $Error=\sqrt{\int |\hat{\mathcal{H}} ~\psi -E~\psi |^2 }$ where the eigen-state energy of equation (2) is calculated as $E=\int \psi^* \hat{\mathcal{H}} ~\psi$ .

The last block ("Schrödinger-equation-component energies") displays the breakdown between the kinetic and potential parts of the ground-state energy.

Plotting the result

At the end of the calculation, the ground-state wave function is stored in the input Schrödinger-equation object SE, together with relevant information such as the domain discretization. For instance, solely relying on SE, one can plot the ground-state wave function with

figure
    plot(SE.xspan,SE.waveFunction.waveFunction,'-','LineWidth',2)
    set(gca,'box','on','FontSize',12,'LineWidth',2)
    xlabel('x (a.u.)')
    ylabel('wave function (a.u.)')
    xlim(SE.xspan([1 end]))

producing (note that the ground-state calculation start from a random seed and thus the resulting wave function is defined with an arbitrary sign that can change from calculation to calculation)

From the plot command line, we see that the domain-discretization grid may be recovered using the xspan property in the object SE (using the standard object-oriented dot notation SE.xspan). On the other hand, the wave function is nested inside another object, which explains the consecutive dots SE.waveFunction.waveFunction. We further explain how to interface input parameters and output results from calculations in the next tutorial.

References

[Javanainen 1988] J. Javanainen, J.H. Eberly, and Q. Su, "Numerical simulations of multiphoton ionization and above-threshold electron spectra," Physical Review A 38, 3430 (1988).

Notes

The results displayed in this tutorial were generated using version 01.21 of the QMol-grid package.

Tutorial design and first-draft writing by William DeNooyer

This wiki is a copy of the documentation provided with the QMol-grid package (accessible in MATLAB documentation, via the "Supplemental Software" section).
Copyright © 2024, Francois Mauger, all right reserved.