QMol_SE - fmauger1/QMol-grid GitHub Wiki

QMol_SE

Schrödinger-equation Hamiltonian and functional.

Description

Use QMol_SE to describe Schrödinger equation (SE) models in simulations. Note that QMol_SE only describes the Schrödinger-equation Hamiltonian and, e.g., does not perform ground state or time-propagation computations. QMol_SE is a handle class.

Class properties

Schrödinger equation model

The QMol_SE class defines the following public get-access properties; each can be changed using the set method:

discretization (disc)

Domain discretization object [ handle object (default []) ]

xspan (x)

Cartesian grid for the domain discretization [ vector (default []) ]

  • Before the QMol_SE object is initialized, as a substitute for the discretization object, xspan is used to directly provide the Cartesian grid (all-increasing, equally-spaced vector) for the domain discretization. Upon initializetion, xspan is converted into the proper discretization object.
  • After the QMol_SE object is initialized, fetches discretization.xspan.
  • Irrespective of the initialization status, xspan consistently provides the Cartesian-grid discretization vector for the SE model.

waveFunction (wfcn)

Wave-function object [ handle object (default []) ]

numberWaveFunction (N)

Number of wave functions in the Schrödinger-equation model [ scalar | vector (default []) ]

  • This is relevant when several wave functions are being considered within the Hamiltonian described by the class.
  • Before the QMol_SE object is initialized, this can be used to set the number of wave functions to include in the Schrödinger-equation model. If left empty, a default of one wave function is assumed.
  • After the QMol_SE object is initialized, fetches the number of wave functions included in the waveFunction component.

potential (V)

Potential object [handle object (default []) ]

Other properties

To facilitate simulations, QMol_SE defines a handful of additional transient properties. These cannot be edited with the set method.

isInitialized (isInit)

Whether the Schrödinger-equation-model object is properly initialized. This is used throughout the QMol-grid package to check that the SE-model object holds meaningful information and is ready for use.

Constant properties

dim (1)

  • Indicates that the Schrödinger-equation model dimension is one.

Class methods

Creation

constructor

Create a Schrödinger-equation-model object with empty class properties.

obj = QMol_SE;

Create a Schrödinger-equation-model object with the name properties set to the specified value. Several name-value pairs can be specified consecutively. Suitable name is any of the Schrödinger-equation-model properties and is case insensitive.

obj = QMol_SE(name1,value1);
obj = QMol_SE(name1,value1,name2,value2,___);

Changing class properties

set

Update the name properties of a Schrödinger-equation-model object to the specified value. Several name-value pairs can be specified consecutively. Suitable name is any of the Schrödinger-equation-model properties and is case insensitive.

obj.set(name1,value1);
obj.set(name1,value1,name2,value2,___);

This is the common name-value pair assignment method used throughout the QMol-grid package. The set method also reset the class. After running, the set property updates the isInitialized flag to a false value.

reset

Reset the object by deleting/re-initializing all run-time properties of the class and updating the isInitialized flag to false.

obj.reset;
  • This is the common reset method available to all classes throughout the QMol-grid package.
  • To avoid any mismatch in the Schrödinger-equation component internal properties, this also resets the discretization (domain discretization), waveFunction (wave functions), and potential (potential) components objects.

clear

Clear all class properties

obj.clear;

Clear a specific set of the class properties. Suitable name is any of the Schrödinger-equation-model properties and is case insensitive.

obj.clear(name1,name2,___);

This is the common clear method available to all classes throughout the QMol-grid package. The clear method also reset the class. The clear method can be used to delete specific properties before saving an instance of the QMol_SE class.

Initializing the object

initialize

Initialize a Schrödinger-equation-model object in silent mode, i.e., without displaying any initialization status, and set the isInitialized flag to true with either

obj.initialize;
obj.initialize(0)
  • This initializes the discretization (domain discretization), waveFunction (wave functions), and potential (potential) component objects.
  • To avoid any mismatch in internal properties, initialize first reset the object before performing the initialization

To perform the same initialization but with displaying the status of the object initialization, including the QMol-grid package header and footer,

obj.initialize(1);

To perform the initialization and display the status of the object initialization, without the QMol-grid package header and footer,

obj.initialize(2);  % 2 or any integer larger than 1

Run-time documentation

getMemoryProfile

Get an estimate of the memory held by a QMol_SE object with either

mem = obj.getMemoryProfile;
mem = obj.getMemoryProfile(false);
  • The object does not need to be initialized for the memory footprint evaluation.
  • getMemoryProfile makes efforts to return the object in the same configuration as before the profiling but small side effects cannot be entirely excluded: all the Schrödinger-equation-model properties are left unchanged while run-time properties (for developers) may be altered. Notably, if the object was not initialized (isInitialized = false ), getMemoryProfile resets the object.
  • The estimate includes (1) the domain discretization, (2) wave functions, and (3) potential. Note that, while the memory estimate tries to be conservative actual (lower or higher) memory requirements may happen in simulations.
  • The output mem is the estimated size in bytes.

Additionally display the detail of the memory footprint with

mem = obj.getMemoryProfile(true);

showDocumentation

Display the run-time documentation for the specific configuration of a QMol_SE object, which must have been initialized beforehand

obj.showDocumentation;

The run-time documentation performs the following steps

  • Display the QMol-grid package header (showing the kernel, external component versions)
  • Display the discretization (domain discretization) object run-time documentation
  • Display the potential (potential) object run-time documentation
  • Display the electronic structure, with the number of wave functions considered
  • Display the bibliography associated with the list of cited references (from the previous steps)
  • Display the funding information
  • Display the QMol-grid package footer

Schrödinger-equation and wave function energies

All the following methods require the object to have been initialized.

getEnergy

Compute the Schrödinger-equation-component energies with either

[Etot,Ekin,Epot] = obj.getEnergy;
[Etot,Ekin,Epot] = obj.getEnergy('Schrodinger equation');
  • Etot is the total (sum of the following components) energy; Ekin is the kinetic energy; Epot is the potential energy. Each of Etot, Ekin, and Epot are the sum of the energy contributions from each of the wave functions in the model.

Compute the wave function energies and error

[E,DE] = obj.getEnergy('wave function');
  • For each wave function $\left|\psi \right\rangle$ , the wave function energy is defined as the real $E=\langle \psi |\hat{\mathcal{H}} |\psi \rangle$ and the error is defined as the positive real $\left\|(\hat{\mathcal{H}} -E)|\psi \rangle \right\|$ . $\hat{\mathcal{H}}$ is the total Schrödinger-equation Hamiltonian operator.
  • This uses the Schrödinger-equation potential.

Specify the potential for computing the wave function energies and error with

[E,DE] = obj.getEnergy('wave function',V);
  • V is a compatible potential

showEnergy

Display the Schrödinger-equation energy components with

obj.showEnergy('Schrodinger equation');

Display the wave function energies and errors with either

obj.showEnergy;
obj.showEnergy('wave function');
  • The energy and error are the same as in getEnergy.

Examples

Create a Schrödinger-equation system

First we define the discretization domain and set of atomic centers in the molecular system

x   =   -20:.1:15;
A1  =   QMol_Va_softCoulomb('name','atom 1','charge',2,'position',-3);
A2  =   QMol_Va_softCoulomb('name','atom 2','charge',2,'position',2);

Next, we define the potential

V   =   QMol_SE_V('atom',{A1,A2});

We can now create the Schrödinger-equation object, initialize it, and access the run-time documentation

SE  =   QMol_SE(                        ...
            'xspan',                x,  ...
            'numberWaveFunction',   3,  ...
            'potential',            V);
SE.initialize;
SE.showDocumentation;

yielding

=== 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 = -20:0.1:15
    Size = 351 (3 x 3 x 3 x 13) 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)
    > atom 1, parameterized as                              (soft Coulomb)
      Z =  2.00 | a =  1.00 | X0 =  -3.00
    > atom 2, parameterized as                              (soft Coulomb)
      Z =  2.00 | a =  1.00 | X0 =   2.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
    3 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).

We can also display the estimated memory footprint for the object

SE.getMemoryProfile(true);
  * Domain discretization                                        
    > domain                                                        2.7 KB
    > gradient                                                      5.5 KB
    > Laplacian                                                     2.7 KB
    > Laplacian (velocity gauge)                                    5.5 KB
  * Wave function(s)                                               16.5 KB
  * Potential                                                    
    > potential                                                     2.7 KB
    > potential gradient                                            2.7 KB

Compute the Schrödinger-equation ground state

To compute the Schrödinger-equation object ground state, we use the QMol_SE_eigs class

GSS = QMol_SE_eigs;
GSS.computeGroundState(SE);

that prints out

=== 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.000 (06/17/2024)                                     F. Mauger

=== Wave-function energies ===============================================
  Wave fcn      Energy (-eV)         Error(a.u.)
  --------     ------------          -----------
      1           52.033              3.473e-13
      2           51.011              3.076e-13
      3           35.681              4.194e-13
  ----------------------------------------------

=== Schrodinger-equation-component energies ==============================
  Component      Energy (a.u.)      Energy (eV)
  -----------    -------------     -------------
  Kinetic             0.703             19.142
  Potential          -5.802           -157.867
  -----------    -------------     -------------
  Total              -5.098           -138.726
  ----------------------------------------------

##########################################################################

To finish with, we plot resulting wave functions

figure; hold on
plot(SE.xspan,SE.wfcn.wfcn(:,1),'-','LineWidth',2,'DisplayName','wfcn #1')
plot(SE.xspan,SE.wfcn.wfcn(:,2),'-','LineWidth',2,'DisplayName','wfcn #2')
plot(SE.xspan,SE.wfcn.wfcn(:,3),'-','LineWidth',2,'DisplayName','wfcn #3')
xlabel('position (a.u.)'); xlim(SE.xspan([1 end]));
ylabel('wave functions')
legend show

Basis-set discretization example

We briefly revisit the Schrödinger-equation model above but using a basis-set discretization. First we define the vectors from which the basis set is to be built

% Domain and atomic centers
x   =   -20:.1:15;
A1  =   QMol_Va_softCoulomb('name','atom 1','charge',2,'position',-3);
A2  =   QMol_Va_softCoulomb('name','atom 2','charge',2,'position',2);

% Common shape for atomic orbitals
AO  =   @(s,n,x0,x) (x(:)-x0).^n .*exp(-(x(:)-x0).^2 * .5/s^2);

% Atomic orbital basis
V   =  [AO(1.3,0,A1.position,x),AO(1.7,1,A1.position,x),AO(0.8,2,A1.position,x), ...
        AO(1.3,0,A2.position,x),AO(1.7,1,A2.position,x),AO(0.8,2,A2.position,x)];

% Orthonormalize the basis
disc=   QMol_disc_basis('xspan',x,'basis',V);
disc.orthonormalizeBasis;

The last command orthonormalizes the atomic-orbital basis, and updates disc with the result.

We can now create the Schrödinger-equation object, like we did in the grid-discretization model above, but specifying the basis-set discretization disc instead of the domain

V   =   QMol_SE_V('atom',{A1,A2});
SE  =   QMol_SE(                            ...
            'discretization',       disc,   ...
            'numberWaveFunction',   3,      ...
            'potential',            V);

Next, we check the model using the run-time documentation

SE.initialize;
SE.showDocumentation;

which gives (we only show the section that is different from above)

(...)
=== Discretization =======================================================
  * Domain discretization                    basis set on a Cartesian grid
    Grid  = -20:0.1:15
    Size  = 351 (3 x 3 x 3 x 13) points
    Basis = 6 vectors
    V-01.21.000 (06/17/2024)                                     F. Mauger
(...)

We can also look at the memory footprint for the basis-set Schrödinger-equation model

SE.getMemoryProfile(true);
  * Domain discretization                                        
    > domain                                                        2.7 KB
    > gradient                                                      5.5 KB
    > Laplacian                                                     2.7 KB
    > Laplacian (velocity gauge)                                    5.5 KB
    > basis                                                        16.5 KB
    > gradient matrix                                               288  B
    > Laplacian matrix                                              288  B
    > Laplacian matrix (velocity gauge)                             576  B
  * Wave function(s)                                                288  B
  * Potential                                                    
    > potential                                                     2.8 KB
    > potential gradient                                            2.8 KB

Note the additional lines for the basis set and associated gradient and Laplacian matrices compared to above. Finally, we compute the ground state using the QMol_SE_eig_basis, which is specifically designed for basis-set Schrödinger-equation models

GSS = QMol_SE_eig_basis;
GSS.computeGroundState(SE);

which yields (once again, we only show the part of the output that are different from above)

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

  * Eigen-state solver for SE Hamiltonians             MATLAB eig function
    using a direct diagonalization of the Hamiltonian matrix.
    V-01.21.000 (06/17/2024)                                     F. Mauger

=== Wave-function energies ===============================================
  Wave fcn      Energy (-eV)         Error(a.u.)
  --------     ------------          -----------
      1           51.977              1.043e-15
      2           50.954              9.504e-16
      3           34.254              7.231e-16
  ----------------------------------------------

=== Schrodinger-equation-component energies ==============================
  Component      Energy (a.u.)      Energy (eV)
  -----------    -------------     -------------
  Kinetic             0.607             16.512
  Potential          -5.648           -153.696
  -----------    -------------     -------------
  Total              -5.041           -137.185
  ----------------------------------------------

##########################################################################

Compared to the grid discretization above, the basis-set wave functions and total Schrödinger-equation energies are all higher (as expected), up to about 1.5 eV. Finally, we plot the new set of orbitals

% Reconstruct orbitals on grid
wfcn =   disc.SE_reconstructWaveFunction(SE.wfcn);

% Plot results
figure; hold on
plot(SE.xspan,wfcn.wfcn(:,1),'-','LineWidth',2,'DisplayName','wfcn #1')
plot(SE.xspan,wfcn.wfcn(:,2),'-','LineWidth',2,'DisplayName','wfcn #2')
plot(SE.xspan,wfcn.wfcn(:,3),'-','LineWidth',2,'DisplayName','wfcn #3')
xlabel('position (a.u.)'); xlim(SE.xspan([1 end]));
ylabel('wave functions')
legend show

(for comparison sake, we slightly edit the plot and add the grid-based result to the plot)

Test suite

For consistency with the rest of the QMol-grid package, QMol_SE defines an associated test suite. Run the test suite for the class in normal or summary mode respectively with

QMol_test.test('SE');
QMol_test.test('-summary','SE');

For developers

For memory profiling (getMemoryProfile), the class also defined the following method (with (Hidden,Access = {?QMol_SE,?QMol_SEq}) attributes)

getDiscCopy

Returns a minimally initialized discretization object

disc = getDiscCopy(obj);
  • If the domain discretization is implicitly defined via x, it is converted into a QMol_disc object, otherwise returns a copy of the Domain discretization. In either cases, the output disc is not initialized to limit memory use.

Notes

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

  • QMol_SE was introduced in version 01.20