20.6. Linear Response Transition Dipole Moment and Related Properties
The current version of PyBEST supports property calculations from pCCD-based (see Detailed input structure of pCCD) Linear Response methods (see The Linear Response Module). Thus, we will focus on a pCCD reference function below. The Transition Dipole Moment (TDM), in the absence of any external potential, can be obtained using only the \(\hat{A}\)-dependent part of the full linear response expression, given as:
Using dipole operators as \(\hat{A}\), one can compute the transition dipole moments for selected excited states.
An intrinsic advantage of using pCCD as the reference wave function, unlike CCSD-based approaches, is that the EOM matrix becomes approximately symmetric. This symmetry allows the transition dipole moment to be approximated using only the right eigenvectors without significantly compromising accuracy.
This approximation enables reduced computational scaling to \(\mathcal{O}(o^2v^2)\).
If you use this module, please cite [ahmadkhani2024].
20.6.1. Implemented Features
Molecular properties are determined from Linear Response calculations. Currently, PyBEST supports two different Linear Response flavors (see The Linear Response Module).
LRpCCD: Linear Response based on a pCCD reference.LRpCCDS: Linear Response based on a pCCD reference with single excitations.
Each class is derived from its corresponding EOM class (see Quick Guide) and reuses most of the functionality.
The corresponding LR property calculations support the following options when determining
transition properties. They are collected as a dictionary and passed as the property_options
argument during function call. Supported dictionary keys are:
- origin:
Coordinates of the reference point for dipole moment (list/tuple containing x, y, z coordinates). Defaults to the center of nuclear charge if not provided.
- transition_dipole_moment:
(bool) Whether to compute transition dipole moments. Default: False.
20.6.2. Performance and Applicability
These LR methods are especially useful for:
Accurate and efficient description of singly and doubly excited states
Cost-effective alternatives to EOM-CCSD for large systems
Systems requiring particle-number-conserving excitation operators
20.6.3. Quick Guide: Transiton dipole moment for LR-pCCD and LR-pCCD+S
This is a short introduction to the Linear Response module. More information on the input and output structure can be found in the following sections. Similar to the previous modules, we will assume the following names for all PyBEST objects
- lf:
a
LinalgFactoryinstance (see Preliminaries).- occ_model:
an Aufbau occupation model of the
AufbauOccModelclass- kin:
the kinetic energy integrals
- ne:
the nucleus-electron attraction integrals
- eri:
the two-electron repulsion integrals
The current version of PyBEST supports two different Linear Response flavors, both based on a
pCCD reference function (The pCCD module). The LRpCCD class allows you to determine
the Linear Response for pair excitations only. A code snippet can be found below.
from pybest.properties import LRpCCD
# For pCCD Linear Response (only electron pairs)
# Compute transition dipole moment
tm_dipole = LRpCCD(lf, occ_model)
results = tm_dipole(
kin,
ne,
eri,
eom_output,
property_options={
"origin": origin_coord,
"transition_dipole_moment": True,
},
printoptions={"nroot": 2},
)
The printoptions dictionary allows you to specify the following option.
- nroot:
(int) Number of excited states to compute (default: 5). Should be the same as in the EOM-pCCD or EOM-pCCD+S calculation.
To include single excitations in the Linear Response, you can use the corresponding
LRpCCDS class.
from pybest.properties import LRpCCDS
# For pCCD+S Linear Response
# Compute transition dipole moment
tm_dipole = LRpCCDS(lf, occ_model)
results = tm_dipole(
kin,
ne,
eri,
eom_output,
property_options={
"origin": origin_coord,
"transition_dipole_moment": True,
},
printoptions={"nroot": 2},
)
The result includes excitation energies and electronic, nuclear, and total transition dipole moments, along with oscillator strength from the EOM matrix.
The results are returned as an IOData container,
while all results are saved to the pybest-results/checkpoint_LR_pCCD.h5
or pybest-results/checkpoint_LR_pCCDS.h5 files.
Specifically, the IOData container contains
the following attributes.
- orb_a:
Molecular orbitals from the pCCD reference
- e_ref:
Total reference energy (pCCD)
- e_ee:
Excitation energies (first entry is 0.0)
- civ_r_ee:
Right eigenvectors (CI-like) of the LR EOM matrix
- t_dm:
Transition dipole moments (if requested)
20.6.4. Example Usage
The example files for obtaining transition dipole moments using the linear response method are collected in the directory data/examples/properties.
20.6.4.1. Example for LR-pCCD transition dipole moments and oscillator strength
This is a basic example of how to calculate properties from linear response methods in PyBEST. This script performs an RHF calculation, followed by a pCCD and EOMpCCD calculation, and an LR-pCCD calculation with a pCCD reference function on the water molecule using the cc-pVDZ basis set.
from pybest import context
from pybest.ee_eom import REOMpCCD
from pybest.gbasis import (
compute_eri,
compute_kinetic,
compute_nuclear,
compute_nuclear_repulsion,
compute_overlap,
get_gobasis,
)
from pybest.geminals import RpCCD
from pybest.linalg import DenseLinalgFactory
from pybest.occ_model import AufbauOccModel
from pybest.properties import LRpCCD
from pybest.wrappers import RHF
# Set up molecule, define basis set
fn_xyz = context.get_fn("test/water.xyz")
basis = get_gobasis("cc-pvdz", fn_xyz)
# Define occupation model, orbitals, and overlap
lf = DenseLinalgFactory(basis.nbasis)
occ_model = AufbauOccModel(basis)
orb_a = lf.create_orbital(basis.nbasis)
olp = compute_overlap(basis)
# Construct Hamiltonian
kin = compute_kinetic(basis)
ne = compute_nuclear(basis)
eri = compute_eri(basis)
nuc = compute_nuclear_repulsion(basis)
nuc = compute_nuclear_repulsion(basis)
# Do a Hartree-Fock calculation
hf = RHF(lf, occ_model)
hf_output = hf(kin, ne, eri, nuc, olp, orb_a)
# Do pCCD calculations
pccd = RpCCD(lf, occ_model)
pccd_output = pccd(kin, ne, eri, hf_output)
# Compute excitation energies using EOM matrix.
eom = REOMpCCD(lf, occ_model)
eom_output = eom(kin, ne, eri, pccd_output, nroot=10)
# Compute transition dipole moments
tm_dipole = LRpCCD(lf, occ_model)
dipole_output = tm_dipole(
kin,
ne,
eri,
eom_output,
property_options={
"origin": (0, 0, 0),
"transition_dipole_moment": True,
},
printoptions={"nroot": 10},
)
20.6.4.2. Example for LR-pCCD+S transition dipole moments and oscillator strength
This is a basic example of how to calculate properties from the LR-pCCD+S method in PyBEST. This script performs an RHF calculation, followed by a pCCD and EOMpCCD+S calculation, and an LR-pCCD+S calculation with a pCCD reference function on the water molecule using the cc-pVDZ basis set.
from pybest import context
from pybest.ee_eom import REOMpCCDS
from pybest.gbasis import (
compute_eri,
compute_kinetic,
compute_nuclear,
compute_nuclear_repulsion,
compute_overlap,
get_gobasis,
)
from pybest.geminals import RpCCD
from pybest.linalg import DenseLinalgFactory
from pybest.occ_model import AufbauOccModel
from pybest.properties import LRpCCDS
from pybest.wrappers import RHF
# Set up molecule, define basis set
fn_xyz = context.get_fn("test/water.xyz")
basis = get_gobasis("cc-pvdz", fn_xyz)
# Define occupation model, orbitals, and overlap
lf = DenseLinalgFactory(basis.nbasis)
occ_model = AufbauOccModel(basis)
orb_a = lf.create_orbital(basis.nbasis)
olp = compute_overlap(basis)
# Construct Hamiltonian
kin = compute_kinetic(basis)
ne = compute_nuclear(basis)
eri = compute_eri(basis)
nuc = compute_nuclear_repulsion(basis)
# Do a Hartree-Fock calculation
hf = RHF(lf, occ_model)
hf_output = hf(kin, ne, eri, nuc, olp, orb_a)
# Do pCCD calculations
pccd = RpCCD(lf, occ_model)
pccd_output = pccd(kin, ne, eri, hf_output)
# Compute excitation energies using EOM matrix.
eom = REOMpCCDS(lf, occ_model)
eom_output = eom(kin, ne, eri, pccd_output, nroot=10)
# Compute transition dipole moments
tm_dipole = LRpCCDS(lf, occ_model)
dipole_output = tm_dipole(
kin,
ne,
eri,
eom_output,
property_options={
"transition_dipole_moment": True,
},
printoptions={"nroot": 10},
)