QMol_CI_conv - fmauger1/QMol-grid GitHub Wiki

`QMol_CI_conv`

CI model, using an explicit convolution scheme to calculate the electron-electron repulsion integrals.

Warning: QMol_CI_conv is a new feature of the QMol-grid package, and has not been fully tested and benchmarked.

Description

Use QMol_CI_conv to describe a configuration-interaction (CI) model associated with a spin-restricted set of orbitals, typically obtained via a density-functional theory (DFT) or Hartree-Fock (HF) ground state calculation. The class also provides functionalities for calculating the CI matrix and its associated ground- and excited-state wave functions, as well as the dipole-coupling matrix.

Class at a glance

Editable properties:

Domain discretization: discretization, xspan
Electronic structure: numberElectron, orbitalBasis, configurationBasis, waveFunction
Potential: externalPotential, interactionPotential
CI calculation: display, tolerance, CImatrix, DXmatrix, isDipole
Configuration basis: type, reference, active, frozen, noDouble, noEmpty

Uneditable properties (read only):

isInitialized, isSetConfigBasis, isSpinPol, dim

Methods:

Changing class properties: set, reset, clear
Initialization: initialize, setConfigurationBasis, cleanConfigurationBasis
Run-time documentation: showDocumentation, getMemoryProfile
CI calculation: computeCImatrix, sparsify, computeGroundState, getEnergy, showEnergy, analyzeCIspectrum
Density: getDensityMatrix, getDensity

Jump to examples

Class properties

CI model

The QMol_CI_conv class defines the following public get-access properties for CI calculation; 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_CI_conv object is initialized, as a substitute for the discretization object, xspan is used to directly provide the Cartesian grid (all-increasing, equally-spaces, vector) for the domain discretization. Upon initialization, xspan is converted into the proper discretization object.
After the QMol_CI_conv object is initialized, fetches discretization``.xspan.
Irrespective of the initialization status, xspan consistently provides the Cartesian-grid discretization vector for the CI model.

`numberElectron (N)`

Number of electrons [ integer| (default []) ]

`orbitalBasis (SOB)`

Basis of spatial orbitals $\phi_1 (r),\phi_2 (r),\ldots,\phi_K (r)$ [ matrix (default []) ]

Each column of orbitalBasis defines the discretization of one spatial orbital basis on the discretization domain.
Spatial orbitals are assumed to form an orthonormal basis $\int \phi_k (r)^* \phi_l (r)~\textrm{d}r=\delta_{kl}$ and no check is performed on this.

`configurationBasis (CSB)`

Configuration-state basis [matrix (default []) ]

Each row of configurationBasis defines a configuration state, as numberElectron vector containing the indexes for the spatial orbitals in orbitalBasis and the sign defining the spin (negative for down and positive for up spins). For instance [-1 -2 1 3] corresponds to the Slater determinant formed by the first two spatial orbitals in orbitalBasis with down spin and orbitals 1 and 3 in orbitalBasis with up spin.
Each row in configurationBasis is assumed to define a unique configuration state and no check is performed on this. Warning: keep in mind that two products of spin orbitals that only differ by the order of the orbitals define linearly dependent configuration states, i.e., [-1 -2 1 3] and [-1 -2 3 1] should not be both included in the configurationBasis; use cleanConfigurationBasis to automatically remove duplicates like these.
In most cases, configurationBasis is built automatically using setConfigurationBasis, to match the Configuration-state basis properties.

`externalPotential (Vext)`

External functional object [ QMol_DFT_Vext object | pseudopotential | cell of pseudopotentials | function handle | vector (default []) ]

A QMol_DFT_Vext object describes the external potential in a similar fashion as for DFT and HF calculations.
A pseudopotential object implicitly describes the external potential as a single pseudopotential.
A cell of pseudopotentials object implicitly describes the external potential as the potential formed by the sum of the provided pseudopotentials.
A function handle specifies the external potential with the signature V = funV(x), where the output V has the same shape as the input vector x and contains element-wise values of the potential at the query points x.
pseudopotential, cell of pseudopotentials, and function handle are converted into a QMol_DFT_Vext object when the CI object is initialized.

`interactionPotential (Vee)`

Electron-electron interaction potential [ function handle (@(x) 1https://github.com/fmauger1/QMol-grid/wiki/sqrt(x.^2+2)) ]

The signature for the (effective) electron-electron interaction potential is V = funV(x), where the output V has the same shape as the input vector x and contains element-wise values of the potential at the query points x.
Because the computations of the electron-electron integrals are performed over an (internally handled) extended domain, a user-defined discretization of the electron-electron interaction potential is not supported (only function handle).

`waveFunction (wfcn)`

Wave function [ [] (default) | column vector | matrix ]

CI wave functions are specified as column vectors of the (complex) amplitudes of the configuration states in the configurationBasis. For example, [1/sqrt(2); 1/sqrt(2); 0; 0; __], with length size(``configurationBasis``,1), has equal population in the first two configuration states.
Several wave functions can be combined in a matrix, which each column corresponding to a different wave function.

`display (disp)`

Whether to display the progress of CI calculations [ true (default) | false ]

display is used by setConfigurationBasis and computeCImatrix.

`tolerance (tol)`

Threshold for setting CI matrix elements to zero [ scalar (default 1e-10) ]

When calculating the CI and dipole-coupling matrices, elements with magnitude smaller than tolerance are forced to zero.
Set a negative tolerance to disable the feature.

`CImatrix (CI)`

CI matrix [ matrix (default []) ]

At the end of a CI calculation, CImatrix contains a copy of the CI matrix.

`isDipole (isDip)`

Whether to calculate the dipole-coupling matrix together with the CI one [ true (default) | false ]

`DXmatrix (DX)`

Dipole-coupling matrix [ matrix (default []) ]

When isDipole is true, at the end of a CI calculation, DXmatrix contains a copy of the dipole coupling matrix.

Configuration-state basis

The QMol_CI_conv class defines the following public get-access properties for automatized building of the configuration-state basis; each can be changed using the set method:

`type`

Type of CI calculation [ 'CIS' (default) | 'CISD' | 'RAS' ]

CIS corresponds to CI with single excitations.
CISD corresponds to CI with single and double excitations. Note: this typically leads to a much bigger configuration state basis.
RAS corresponds to CI with a restricted active space. If no noDouble or noEmpty are specified, this corresponds to a complete active space (CAS). Warning: the size of the configuration state basis grows exponentially fast with the size of the active space.
type is used by setConfigurationBasis when building configurationBasis.

`reference (ref)`

Configuration state reference [ vector (default []) ]

numberElectron vector containing the indexes for the spatial orbitals in orbitalBasis and the sign defining the spin (negative for down and positive for up spins). For instance [-1 -2 1 3] corresponds to the Slater determinant formed by the first two spatial orbitals in orbitalBasis with down spin and orbitals 1 and 3 in orbitalBasis with up spin.
If left empty, a spin-restricted system is assumed and the first numberElectron``/2 spatial orbitals in orbitalBasis are used for the reference. If numberElectron is an odd number an error is throwed when building the configurationBasis.
reference is used by setConfigurationBasis when building configurationBasis.
For CIS(D) calculations, reference defines the spin orbitals from which the excitations should be considered.
For RAS calculations, a reference is still required and used to define the number of electrons to include in each spin channel. For instance [-1 -2 1 3 4] indicates that 2 electrons should be put in the down-spin channel and 3 in the up one.

`active (act)`

List of spin orbitals that may contain an electron [ vector (default []) ]

Vector containing the "active" spin orbitals in the CIS(D)/RAS basis, defined as the indexes for the spatial orbitals in orbitalBasis and the sign indicating the spin (negative for down and positive for up spins). For instance [-1 -2 -3 -4 1 2 3] indicates that the spin orbitals formed by the first 4 spatial orbitals with down spin and only the first 3 spatial orbitals with up spin should be considered for building the configurationBasis.
The largest index magnitude in the active vector must not be larger than the number of spatial orbitals in orbitalBasis.
If left empty, all the spatial orbitals in orbitalBasis with both up and down spins are included in the active space.
active is used by setConfigurationBasis when building configurationBasis.
For CIS(D) calculations, only active spin orbitals that are different from the reference are relevant for building building the configurationBasis. Duplicates with the reference are ignored (including the in the run-time documentation).
For RAS calculations, all active spin orbitals must be listed as the reference is only used to determine the number of electrons in each spin channel.

`frozen (frz)`

List of spin orbitals that must contain an electron [ vector (default []) ]

Vector containing the "frozen" spin orbitals, which are included in all the vectors in the CIS(D)/RAS basis. Frozen orbitals are defined via the indexes for the spatial orbitals in orbitalBasis and the sign indicating the spin (negative for down and positive for up spins). For instance [-1 1 2] indicates all configuration states in the configurationBasis must include the spin orbitals formed by the first spatial orbitals with down spin and the first 2 spatial orbitals with up spin.
The largest index magnitude in the frozen vector must not be larger than the number of spatial orbitals in orbitalBasis.
If left empty, no orbitals are frozen.
frozen is used by setConfigurationBasis when building configurationBasis.
For CIS(D) calculations, all frozen orbitals should also be included in the reference.
For RAS calculations, the number of frozen orbitals in each spin channel must be compatible with the number of electron each hold (defined via the reference).

`noDouble (noDbl)`

List of spatial orbitals that must not contain two electrons [ vector (default) ]

Vector containing the indexes of the spatial orbitals that may be empty or contain one electron but should not hold two electrons. For instance [5 6] indicates that spatial orbitals 5 and 6 may not be fully occupied.
If left empty, the feature is disabled.
noDouble is used by setConfigurationBasis when building configurationBasis.
For CIS(D) calculations, the constraint is not imposed on the reference configuration.
After building the CIS(D)/RAS basis associated with the reference, active, and frozen spin orbitals, configuration states that allocate two electrons in the noDouble spatial orbitals are removed from the configurationBasis.

`noEmpty (noEmp)`

List of spatial orbitals that must not be empty [ vector (default) ]

Vector containing the indexes of the spatial orbitals that may contain one or two electrons but should not be empty. For instance [1 2] indicates that the first two spatial orbitals should always hold at least one electron.
If left empty, the feature is disabled.
noEmpty is used by setConfigurationBasis when building configurationBasis.
For CIS(D) calculations, the constraint is not imposed on the reference configuration.
After building the CIS(D)/RAS basis associated with the reference, active, and frozen spin orbitals, configuration states that do not allocate any electron in the noEmpty spatial orbitals are removed from the configurationBasis.

Other properties

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

`isInitialized (isInit)`

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

`isSetConfigBasis`

Tacks whether the configuration basis configurationBasis was build using [what] or user defined [ true | false (default) ]

Constant properties

`isSpinPol (true)`

Indicates that, by construction, the CI model is spin polarized.

`dim (1)`

Indicates that the CI model dimension is one.

Class methods

Creation

constructor

Create a CI-model object with empty class properties.

obj = QMol_CI_conv;

Create a CI-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 CI-model and Configuration-state basis properties, and is case insensitive.

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

Create a CI-model object with the discretization (domain discretization), numberElectron (number of electrons), orbitalBasis (spatial-orbital basis), externalPotential (external potential), and interactionPotential (electron-electron interaction potential) properties are imported from a QMol_DFT_spinRes DFT object by passing it as the first argument of the constructor. Additional or overwriting parameters can be specified as a list of name-value pair arguments.

obj = QMol_CI_conv(DFT);
obj = QMol_CI_conv(DFT,name1,value1,name2,value2,___);

For basis-set DFT and HF models, orbitalBasis contains the reconstructed spatial orbitals over the discretization domain.

Changing class properties

`set`

Update the name properties of a CI-model object to the specified value. Several name-value pairs can be specified consecutively. Suitable name is any of the CI-model and Configuration-state basis properties, and is case insensitive.

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

Update the discretization (domain discretization), numberElectron (number of electrons), orbitalBasis (spatial-orbital basis), externalPotential (external potential), and interactionPotential (electron-electron interaction potential) properties from a QMol_DFT_spinRes DFT object by passing it as the first argument. Additional or overwriting parameters can be specified as a list of name-value pair arguments.

obj.set(DFT);
obj.set(DFT,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 and isSetConfigBasis flags 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 CI component internal properties, this also resets the discretization (domain discretization) and externalPotential (external potential) components objects.

`clear`

Clear all class properties.

obj.clear;

Clear a specific set of the class properties. Suitable name is any of the CI-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_CI_conv class.

Initializing the object

`initialize`

Initialize a CI-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) and externalPotential (external potential)component objects as well as the internal discretization of interactionPotential on the extended domain used to calculate electron-electron interaction integrals.
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 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

`setConfigurationBasis`

Set configurationBasis to match the Configuration-state basis properties.

obj.setConfigurationBasis;

setConfigurationBasis is typically used between the creation of the CI object and the calculation of the CI matrix, using computeCImatrix.
setConfigurationBasis first (re)initializes the CI-model object before building configurationBasis.
setConfigurationBasis sets isSetConfigBasis to true, to indicate that configurationBasis matches the Configuration-state basis properties.
If display is true, setConfigurationBasis also prints out a summary of the configuration state basis.

`cleanConfigurationBasis`

Clean the configuration state basis by removing duplicate configuration states in configurationBasis, including those related by a permutation of the spin orbitals in the index vector.

obj.cleanConfigurationBasis;

The cleaning process preserves the order of (i) configuration states in the basis, keeping only the first unique instance of each state, and (ii) the order of the spin-orbital withing each configuration state.
For instance, the configurationBasis [-1 -2 1 2; -1 -3 1 2; -2 -1 1 2;-1 -2 1 3] is cleaned to [-1 -2 1 2; -1 -3 1 2; -1 -2 1 3] since [-2 -1 1 2] is a duplicate of [-1 -2 1 2], and thus removed from the set.

Run-time documentation

`getMemoryProfile`

Get an estimate of the memory held by a QMol_CI_conv 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 CI-model and Configuration-state basis properties are left unchanged while run-time properties (for developers) may be altered.
The estimate includes the (1) domain discretization, (2) external potential, (3) electron-electron interaction potential, (4) spatial orbital basis, (5) configuration state basis, (6) CI matrix, and (7) dipole-coupling matrix if isDipole is set to true. Note that all these components may not be used in actual simulations and the memory estimate tries to be conservative. On the other hand, it only includes the discretization of member components on the domain grid and ignores other (small) properties.
If configurationBasis is empty, the profiling assumes that the configuration state basis will be built using the setConfigurationBasis method. Otherwise, the profiling uses configurationBasis to estimate the memory footprint.
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_CI_conv 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 electronic structure
Display the configuration state basis properties
Display the externalPotential (external functional and potential) object run-time documentation
Display the HartreePotential object run-time documentation
Display the interactionPotential
Display the bibliography associated with the list of cited references (from the previous steps)
Display the funding information
Display the QMol-grid package footer

CI calculation

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

`computeCImatrix`

Computes the CI matrix associted with the CI-model properties.

obj.computeCImatrix;

The calculation first precomputes the core Hamiltonian matrix and two-electron 4-tensor, which are then used to build the CI matrix.
The QMol_CI_conv object must have a properly defined configurationBasis, e.g., built with setConfigurationBasis, to perform the CI matrix calculation.
The CI matrix is stored in the CImatrix member property. If isDipole is true, the dipole-coupling matrix is also calculated and stored in the DXmatrix member property. Both CImatrix and DXmatrix are size(``configurationBasis``,1)-by-size(``configurationBasis``,1) matrices.
If display is true, computeCImatrix also prints out the progress of the CI-matrix calculation.
In the core Hamiltonian, CI, and dipole-coupling matrices as well as the two-electron 4-tensor, elements with magnitude smaller than tolerance are forced to zero (set tolerance to a negative value to disable this feature).

Get a copy of the CI matrix with

CI = obj.computeCImatrix;

When isDipole is true, get a copy of the dipole-coupling matrix with

[CI,DX] = obj.computeCImatrix;

`sparsify`

Convert the CI and dipole-coupling into sparse matrices

obj.sparsify;

This can accelerate some calculations when the CImatrix is large.

`computeGroundState`

Compute the CImatrix ground state with

obj.computeGroundState;
obj.computeGroundState(1);

The result is stored in waveFunction.
The CI object must contain a valid CImatrix.

Additionally compute a set of lowest-energy excited states with

obj.computeGroundState(n);

The integer n specifies the number of states to calculate, thus corresponding to n-1 excited states.

`getEnergy`

Compute the member wave function energy and error

[E,DE] = obj.getEnergy;

For each wave function $|\psi \rangle$ in waveFunction, the energy is defined as the real $E=\langle \psi |\hat{\mathcal{H}} |\psi \rangle$ and the error as the positive real $\|(\hat{\mathcal{H}} -E)|\psi \rangle \|$ , where $\hat{\mathcal{H}}$ is the CI Hamiltonian operator (defined via the CImatrix).
The CI object must contain a valid CImatrix.

`showEnergy`

Display the member wave function energy and error

obj.showEnergy;

The energy and error are obtained from getEnergy.

`analyzeCIspectrum`

After calculating the CI matrix, use analyzeCIspectrum to print out the energies and wave function composition of eigen states of the CI model. Display the analysis for the n lowest energy eigen states with

obj.analyzeCIspectrum(n);

For each eigen state, analyzeCIspectrum prints out the total electronic energy for the state as well as the configuration states in configurationBasis that substantially contribute to the state (with population larger than 1%).
For excited state, analyzeCIspectrum also prints out the excitation energy, defined as the energy difference with the ground state.

One-body density and one-particle density matrices

`getDensity`

Compute the one-body density associated with a wave function of a CI model, defined as $\rho (r)=\int |\psi (r,\omega ,x_2 \ldots x_N )|^2 ~d\omega dx_2 \ldots dx_N$ , with the integration running over (i) the first electron spin coordinate and (ii) all the following electronic and spin coordinates. The one-body density is further decomposed between its up- and down-spin contributions.

rho = obj.getDensity;

Returns the one-body density of the waveFunction property, which must contain a single wave function.
The output rho is a QMol_DFT_density object containing the discretization of the up- and down-spin one-body density components. The density object is initialized.
Internally, the one-particle reduced density matrices associated with psi are calculated, using getDensityMatrix, from which the one-body density is reconstructed.

rho = obj.getDensity(psi);

Returns the one-body density for the input wave function psi instead of waveFunction.
The input wave function psi is specified as a column vector of population coefficients in the configuration state basis configurationBasis.

Instead of the wave function, the one-body density can be defined via the one-particle density matrices

rho = obj.getDensity(rhoUp,rhoDown);

The two input rhoUp and rhoDown matrices are typically obtained using getDensityMatrix.

Optionally specify the QMol_DFT_density object in which to store the one-body density with

obj.getDensity(___,rho);

`getDensityMatrix`

Computes the one-particle reduced density matrix (1-RDM) associated with the wave function(s) of a CI model. The matrix elements $\rho_{m,n}$ of the 1-RDM are associated with the one-particle density operators $|m\rangle \langle n|$ of the $m^{\textrm{th}}$ and $n^{\textrm{th}}$ spin orbitals in the basis. The one-body density $\rho$ is recovered from the 1-RDM as $\rho =\sum_{m,n} \rho_{m,n} \chi_m^* \chi_n$ . The 1-RDM is further decomposed between its up- and down-spin contributions (since the two spin channels do not mix in the joint 1-RDM space).

[rhoUp,rhoDown] = obj.getDensityMatrix;

Returns the 1-RDM matrices for the up (rhoUp) and down (rhoDown) spin channels of the waveFunction property.
For a single wave function (column vector waveFunction), the 1-RDM are size(``configurationBasis``,1)-by-size(``configurationBasis``,1) matrices.
If waveFunction contains several wave functions, the output rhoUp and rhoDown are size( configurationBasis``,1)-by-size(``configurationBasis``,1)- by-size(``waveFunction``,2) arrays. One can recover the 1-RDM matrices for the k $^{\textrm{th}}$ wave function component with rhoUp(:,:,k) and rhoDown(:,:,k), which are all Hermitian matrices.
To speed up calculations of the density matrices, contributions smaller than tolerance are skipped in 1-RDM calculation. This may lead to slight difference when the 1-RDM for a collection of wave functions are calculated together (matrix psi) vs. separately (vector psi). Set tolerance to a negative value to disable this feature.
Note: the calculation of the 1-RDM can be time consuming, especially when using large configuration state bases (the algorithm shares a similar outline as when calculating the CI matrix).

[rhoUp,rhoDown] = obj.getDensityMatrix(psi);

Returns the 1-RDM matrices for the input wave function(s) psi instead of waveFunction.
Each column in psi specifies the population coefficients in the configuration state basis configurationBasis.

Examples

Build the spatial orbital basis set

To define the spatial-orbital basis set for our CI calculations, we first run a spin-restricted DFT calculation

% Molecular model
x       =   -50:.1:50;
A1      =   QMol_Va_softCoulomb('name','atom 1','charge',3,'position',-3);
A2      =   QMol_Va_softCoulomb('name','atom 2','charge',3,'position',2);

% DFT model
Vext    =   QMol_DFT_Vext('atom',{A1,A2});
Vh      =   QMol_DFT_Vh_conv;
Vxc     =   {QMol_DFT_Vx_LDA_soft,QMol_DFT_Vc_LDA_soft};

DFT     =   QMol_DFT_spinRes(                                   ...
                'xspan',                        x,              ...
                'occupation',                   [2 2 2 0 0 0],  ...
                'externalPotential',            Vext,           ...
                'HartreePotential',             Vh,             ...
                'exchangeCorrelationPotential', Vxc,            ...
                'selfInteractionCorrection',    'ADSIC'         );

Note that we include 3 virtual orbitals (with occupation 0) in the DFT model, which will define the active space for the CI calculations. Finally, we run a DFT-ground state calculation to obtain our spatial orbital basis set.

SCF     =   QMol_DFT_SCF_Anderson;
SCF.solveSCF(DFT);

yielding

(...)
=== Orbital energies =====================================================
  Orbital      Occ. (elec.)         Energy (-eV)               Error(a.u.)
  -------      ------------         ------------               -----------
      1            2.00                29.364                   2.051e-12
      2            2.00                28.874                   2.104e-12
      3            2.00                10.650                   4.065e-12
      4            0.00                 7.672                   8.120e-12
      5            0.00                 3.807                   1.381e-11
      6            0.00                 2.101                   1.294e-11
  ------------------------------------------------------------------------

=== DFT-component energies ===============================================
  Component      Energy (a.u.)      Energy (eV)
  -----------    -------------     -------------
  Kinetic             1.352             36.802
  External          -17.135           -466.274
  Hartree             5.521            150.243
  Exch.-corr.        -0.557            -15.170
  -----------    -------------     -------------
  Total             -10.819           -294.398
  ----------------------------------------------

Set the CI model

We now aim to define our CI model associated with the basis of DFT spatial orbitals above. To do so, we let QMol_CI_conv read the relevant properties off the DFT object.

CI      =   QMol_CI_conv(DFT);

Next, we aim to perform a CISD calculation using the occupied DFT orbitals as our reference and all the virtual DFT orbitals for the active space. These correspond to the default values for the reference and active properties, and thus only need to set the type of CI.

CI.set('type','CISD');

Before building the configuration space basis, we check the memory footprint of our CI model with

CI.getMemoryProfile(true);

that prints out

  * Domain discretization                                        
    > domain                                                        7.8 KB
    > gradient                                                     15.6 KB
    > Laplacian                                                     7.8 KB
    > Laplacian (velocity gauge)                                   15.6 KB
  * External functional                                          
    > potential                                                     7.8 KB
    > potential gradient                                            7.8 KB
  * Electron-electron interaction                                  15.6 KB
  * Spatial orbital basis                                          46.9 KB
  * Configuration state basis                                       5.5 KB
  * CI matrix                                                     108.8 KB
  * Dipole-coupling matrix                                        108.8 KB

Finally, we build the CISD basis with

CI.setConfigurationBasis;

that prints out

(...)
=== Build the configuration-interaction (CI) model [Visentin 2025] =======
  * Discretization                                                      OK
  * External (atomic or molecular) potential                            OK
  * Electron-electron interaction potential                             OK
=== References ===========================================================
  [Mauger 2024b] F. Mauger and C. Chandre, "QMol-grid : A MATLAB package
    for quantum-mechanical simulations in atomic and molecular systems," 
    SoftwareX 28, 101968 (2024).
  [Visentin 2025] G. Visentin et al., in preparation (2025).

=== Funding ==============================================================
    The original development of the QMol-grid package, 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 and CI features was supported by the National
  Science Foundation under Grant No. PHY-2207656.

##########################################################################
=== Set the configuration-state basis ====================================
  * Initialization                                                    DONE
  * Configuration-state basis
    Type         = CISD
    Total spin   = 0
    Reference(s) = -1 -2 -3 1 2 3
    Active       = 4 5 6
                   -6 -5 -4
    118 configuration-state basis elements
    V-01.22.000 (05/20/2025)                       G. Visentin & F. Mauger

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

Compute the CI matrix and analyze its spectrum

We now have all the elements in place to calculate the CI matrix with

CI.computeCImatrix;

that prints out

(...)
=== Build the configuration-interaction (CI) model [Visentin 2025] =======
  * Discretization                                                      OK
  * External (atomic or molecular) potential                            OK
  * Electron-electron interaction potential                             OK

=== References ===========================================================
  [Mauger 2024b] F. Mauger and C. Chandre, "QMol-grid : A MATLAB package
    for quantum-mechanical simulations in atomic and molecular systems," 
    SoftwareX 28, 101968 (2024).
  [Visentin 2025] G. Visentin et al., in preparation (2025).

=== Funding ==============================================================
    The original development of the QMol-grid package, 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 and CI features was supported by the National
  Science Foundation under Grant No. PHY-2207656.

##########################################################################
=== CI-matrix calculation ================================================
                                       Calculation progress (%)
  Calculation step          0       20       40       60       80      100
  ----------------------    ----------------------------------------------
  Core integrals            ||||||||||||||||||||||||||||||||||||||||||||||
  Dipole coupling elements  ||||||||||||||||||||||||||||||||||||||||||||||
  Two-electron integrals    ||||||||||||||||||||||||||||||||||||||||||||||
  CI matrix                 ||||||||||||||||||||||||||||||||||||||||||||||
  ------------------------------------------------------------------------
                                            CI-matrix calculation complete

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

Finally, we inspect the ground and first 2 excited states of our CI model with

CI.analyzeCIspectrum(3);

that prints out

=== CI matrix spectrum ===================================================
  * CI matrix
    size = 118 x 118 with 3374 nonzero elements
  * Ground state
    Total energy      =    -10.838 a.u. =   -294.920 eV
    >  92.581 % [  -1  -2  -3  1  2  3  ]
    >   6.250 % [  -1  -2  -4  1  2  4  ]
  * Excited state
    Total energy      =    -10.769 a.u. =   -293.039 eV
    Excitation energy =      0.069 a.u. =      1.881 eV
    >  49.018 % [  -1  -2  -4  1  2  3  ]
    >  49.018 % [  -1  -2  -3  1  2  4  ]
  * Excited state
    Total energy      =    -10.623 a.u. =   -289.062 eV
    Excitation energy =      0.215 a.u. =      5.857 eV
    >   3.828 % [  -1  -2  -6  1  2  3  ]
    >  38.990 % [  -1  -2  -4  1  2  3  ]
    >  38.990 % [  -1  -2  -3  1  2  4  ]
    >   3.828 % [  -1  -2  -3  1  2  6  ]
    >   7.014 % [  -1  -2  -5  1  2  4  ]
    >   7.014 % [  -1  -2  -4  1  2  5  ]
==========================================================================

In the ground state, we notice a substantial (>6%) contribution from the double excitation from the highest-occupied (HOMO, index +/- 3) to the lowest unoccupied (LUMO, index +/-4) molecular orbitals. The first excited state sits less than 2 eV above the ground one, and is almost fully composed of single excitations from the HOMO to the LUMO orbitals (one in each spin channel). On the other hand, the second excited state contains substantial contributions from single excitations from the HOMO to LUMO and LUMO+2 (index +/- 6), as well as double excitations from the HOMO to the LUMO and LUMO+1 (index +/-5).

Plot the ground state one-body density

To get a better visualization of the ground state, we now want to plot its one-body density. First, we calculate the ground state vector, which corresponds to the eigen vector of the CI matrix with the smallest energy

CI.computeGroundState;

from which we can recover the one-body density (and that of the parent DFT model for comparison)

rho_CI  = CI.getDensity(CI.waveFunction);
rho_DFT = DFT.getDensity;

Finally, we plot the results

figure
    subplot(3,1,[1 2]), hold on
        plot(CI.xspan,rho_CI.densityUp+rho_CI.densityDown,'-','LineWidth',2,'DisplayName','CI')
        plot(DFT.xspan,rho_DFT.density,':','LineWidth',2,'DisplayName','DFT')
        xlim(DFT.xspan([1 end]));
        ylabel('densty (elechttps://github.com/fmauger1/QMol-grid/wiki/a.u.)')
        legend show
    subplot(3,1,3), hold on
        plot(CI.xspan,rho_CI.densityUp+rho_CI.densityDown-rho_DFT.density,'-','LineWidth',2)
        xlabel('position (a.u.)'); xlim(DFT.xspan([1 end]));
        ylabel('difference')

Basis-set discretization example

To illustrate CI calculations starting from a basis-set discretization, we consider a HF set of orbitals obtained with

% Domain and atomic centers
x       =   -20:.1:15;
A1      =   QMol_Va_softCoulomb('name','atom 1','charge',3,'position',-3);
A2      =   QMol_Va_softCoulomb('name','atom 2','charge',3,'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;

% HF model
Vext    =   QMol_DFT_Vext('atom',{A1,A2});
Vh      =   QMol_DFT_Vh_conv;
Vxc     =   QMol_DFT_Vx_XX_conv;

HF      =   QMol_DFT_spinRes(                                   ...
                            'discretization',               disc,           ...
                            'occupation',                   [2 2 2 0 0 0],  ...
                            'externalPotential',            Vext,           ...
                            'HartreePotential',             Vh,             ...
                            'exchangeCorrelationPotential', Vxc);

% DFT ground state
SCF     =   QMol_DFT_SCF_Anderson;
SCF.solveSCF(HF);

yielding

(...)
=== Orbital energies =====================================================
  Orbital      Occ. (elec.)         Energy (-eV)               Error(a.u.)
  -------      ------------         ------------               -----------
      1            2.00                34.843                   3.274e-12
      2            2.00                34.125                   2.177e-12
      3            2.00                 9.103                   4.466e-12
      4            0.00                -0.350                   3.637e-12
      5            0.00               -48.715                   3.166e-12
      6            0.00               -49.452                   3.523e-12
  ------------------------------------------------------------------------

=== DFT-component energies ===============================================
  Component      Energy (a.u.)      Energy (eV)
  -----------    -------------     -------------
  Kinetic             1.262             34.338
  External          -16.974           -461.887
  Hartree             6.639            180.650
  Exch.-corr.        -1.652            -44.945
  -----------    -------------     -------------
  Total             -10.725           -291.844
  ----------------------------------------------

Note that the LUMO+1 and LUMO+2 orbitals have very large energies and may not make much sense to include in practical CI calculations. We keep them here for illustration purposes.

Next, we set the CI model with a complete active space (RAS type with no additional restrictions) over all the calculated HF molecular orbitals and check the memory footprint of the system

CI      =   QMol_CI_conv(HF,'type','RAS');
CI.getMemoryProfile(true);

yielding

(...)
  * 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
  * External functional                                          
    > potential                                                     2.7 KB
    > potential gradient                                            2.7 KB
  * Electron-electron interaction                                   5.5 KB
  * Spatial orbital basis                                          16.5 KB
  * Configuration state basis                                      18.8 KB
  * CI matrix                                                       1.2 MB
  * Dipole-coupling matrix                                          1.2 MB

Note that the HF basis-set discretization domain has been converted to a regular discretization grid as the spatial-orbital basis is rebuilt over the discretization grid in QMol_CI_conv.

From there, the steps are identical to the grid-based example above. We first build the configuration state basis with

CI.setConfigurationBasis;

that prints out

(...)
=== Build the configuration-interaction (CI) model [Visentin 2025] =======
  * Discretization                                                      OK
  * External (atomic or molecular) potential                            OK
  * Electron-electron interaction potential                             OK
=== References ===========================================================
  [Mauger 2024b] F. Mauger and C. Chandre, "QMol-grid : A MATLAB package
    for quantum-mechanical simulations in atomic and molecular systems," 
    SoftwareX 28, 101968 (2024).
  [Visentin 2025] G. Visentin et al., in preparation (2025).

=== Funding ==============================================================
    The original development of the QMol-grid package, 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 and CI features was supported by the National
  Science Foundation under Grant No. PHY-2207656.

##########################################################################
=== Set the configuration-state basis ====================================
  * Initialization                                                    DONE
  * Configuration-state basis
    Type         = RAS
    Total spin   = 0
    Active       = 1 2 3 4 5 6
                   -6 -5 -4 -3 -2 -1
    400 configuration-state basis elements
    V-01.22.000 (05/20/2025)                       G. Visentin & F. Mauger

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

Next, we calculate the CI matrix with

CI.computeCImatrix;

that prints out

(...)
=== Build the configuration-interaction (CI) model [Visentin 2025] =======
  * Discretization                                                      OK
  * External (atomic or molecular) potential                            OK
  * Electron-electron interaction potential                             OK

=== References ===========================================================
  [Mauger 2024b] F. Mauger and C. Chandre, "QMol-grid : A MATLAB package
    for quantum-mechanical simulations in atomic and molecular systems," 
    SoftwareX 28, 101968 (2024).
  [Visentin 2025] G. Visentin et al., in preparation (2025).

=== Funding ==============================================================
    The original development of the QMol-grid package, 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 and CI features was supported by the National
  Science Foundation under Grant No. PHY-2207656.

##########################################################################
=== CI-matrix calculation ================================================
                                       Calculation progress (%)
  Calculation step          0       20       40       60       80      100
  ----------------------    ----------------------------------------------
  Core integrals            ||||||||||||||||||||||||||||||||||||||||||||||
  Dipole coupling elements  ||||||||||||||||||||||||||||||||||||||||||||||
  Two-electron integrals    ||||||||||||||||||||||||||||||||||||||||||||||
  CI matrix                 ||||||||||||||||||||||||||||||||||||||||||||||
  ------------------------------------------------------------------------
                                            CI-matrix calculation complete

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

Finally, we inspect the ground and first 2 excited states of our CI model with

CI.analyzeCIspectrum(3);

that prints out

=== CI matrix spectrum ===================================================
  * CI matrix
    size = 400 x 400 with 24432 nonzero elements
  * Ground state
    Total energy      =    -10.758 a.u. =   -292.741 eV
    >  91.904 % [  -3  -2  -1  1  2  3  ]
    >   7.859 % [  -4  -2  -1  1  2  4  ]
  * Excited state
    Total energy      =    -10.698 a.u. =   -291.114 eV
    Excitation energy =      0.060 a.u. =      1.627 eV
    >  49.753 % [  -4  -2  -1  1  2  3  ]
    >  49.753 % [  -3  -2  -1  1  2  4  ]
  * Excited state
    Total energy      =    -10.499 a.u. =   -285.686 eV
    Excitation energy =      0.259 a.u. =      7.055 eV
    >  49.766 % [  -4  -2  -1  1  2  3  ]
    >  49.766 % [  -3  -2  -1  1  2  4  ]
==========================================================================

Without much surprise, we see that none of the ground and first 3 excited state involve the LUMO+1 (index +/-5) and LUMO+2 (index +/- 6), which are very far apart energetically (see the comment above).

Test suite

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

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

For developers

For internal use QMol_CI_conv defines additional transient properties. These cannot be edited with the set method but can be directly edited by classes in the QMol-grig package.

V2e

Internal electron-electron interaction potential [ vector (default []) ]

Discretization of interactionPotential over the extended domain used to calculate electron-electron interaction integrals.
This property is set when the CI object is initialized.

`isReal`

Whether the spatial orbitals are real valued [ true | false | (default []) ]

This property is set when the CI object is initialized.

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 case, the output disc is not initialized to limit memory use.

getConfigurationBasisParameters

Parsing of the flavor of CI (CIS, CISD, or RAS) and their associated number of electrons, reference, active space, and frozen spin orbitals. It provides a common interface for building the configuration-space basis and run-time documentation.

[algo,locN,locRef,locAct,locFrz] = obj.getConfigurationBasisParameters(isBld)

algo = 1 corresponds to CIS, algo = 2 corresponds to CISD, and algo = 0 corresponds to RAS.
locN is the number of electrons.
locRef are the indexes for the spin orbitals forming the reference configuration state.
locAct are the indexes for the spin orbitals forming the active space.
locFrz are the indexes for the frozen spin orbitals.
Set isBld = true when called inside setConfigurationBasis and isBld = false otherwise for proper run-time documentation display.

Notes

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

QMol_CI_conv was introduced in version 01.22