.. : PyBEST: Pythonic Black-box Electronic Structure Tool : Copyright (C) 2016-- The PyBEST Development Team : : This file is part of PyBEST. : : PyBEST is free software; you can redistribute it and/or : modify it under the terms of the GNU General Public License : as published by the Free Software Foundation; either version 3 : of the License, or (at your option) any later version. : : PyBEST is distributed in the hope that it will be useful, : but WITHOUT ANY WARRANTY; without even the implied warranty of : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the : GNU General Public License for more details. : : You should have received a copy of the GNU General Public License : along with this program; if not, see : -- .. _user_properties_lr: Linear Response Transition Dipole Moment and Related Properties ############################################################### The current version of PyBEST supports property calculations from pCCD-based (see :ref:`pccd`) Linear Response methods (see :ref:`user_lr_intro`). 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 :math:`\hat{A}`-dependent part of the full linear response expression, given as: .. math:: \Gamma_{0 \rightarrow k}^{\hat{A}} &= \langle 0 | \hat{A} | k \rangle \\ &= \langle \Lambda |[\hat{A}, \, \hat{\tau}_k]|\mathrm{pCCD}\rangle \\ &\quad - \sum_{\nu} (-\hat{J} + \omega \boldsymbol{I})_{\mu\nu}^{-1} \langle \bar{\nu} | \hat{A} | \mathrm{pCCD} \rangle \langle \Lambda | [[\hat{H}_0, \hat{\tau}_\nu], \hat{\tau}_k] | \mathrm{pCCD} \rangle Using dipole operators as :math:`\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 Jacobian 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 :math:`\mathcal{O}(o^2v^2)`. If you use this module, please cite [ahmadkhani2024]_. Implemented Features ==================== Molecular properties are determined from Linear Response calculations. Currently, PyBEST supports two different Linear Response flavors (see :ref:`user_lr_intro`). - :class:`~pybest.properties.lr_pccd.LRpCCD`: Linear Response based on a pCCD reference - :class:`~pybest.properties.lr_pccd_s.LRpCCDS`: Linear Response based on a pCCD with single excitations Each class is derived from its corresponding EOM class (see :ref:`user_eomcc_intro`) and reuses most of the functionality, modifying only the construction of the Hamiltonian to match the true LR structure. 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: :operator_A: One-electron operator (e.g., dipole integrals) :operator_B: One-electron operator (e.g., dipole integrals) :coordinates: Used for origin-dependent properties (list containing x, y, z coordinates) :transition_dipole_moment: (bool) Whether to compute transition dipole moments 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 Quick Guide: Transiton dipole moment for pCCD and 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 :py:class:`~pybest.linalg.base.LinalgFactory` instance (see :ref:`user_linalg_intro`). :occ_model: an Aufbau occupation model of the :py:class:`~pybest.scf.occ.AufbauOccModel` class :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 (:ref:`user_pccd`). For the transition dipole moment calculation, you need to pass the integrals of the dipole operator (see :ref:`user_molecularham_dipole`) and the molecular coordinates of the origin. The code snippet below demonstrates how to generate or extract this data from previously generated PyBEST output ``pccd_output``. .. code-block:: python from pybest.gbasis import compute_dipole from pybest.utility import get_com # Compute dipole integrals for center of charge/mass x, y, z = get_com(basis) dipole = compute_dipole(basis, x=x, y=y, z=z) # Return the origin from previously generated output from pybest.wrappers.multipole import check_coord # check_coord returns the origin of the dipole integrals coord = check_coord(dipole, pccd_output) The :py:class:`~pybest.properties.lr_pccd.LRpCCD` class allows you to determine the Linear Response for pair excitations only. A code snippet can be found below. .. code-block:: python from pybest.properties import LRpCCDS # For pCCD Linear Response (only electron pairs) # Compute transition dipole moment tm_dipole = LRpCCD(lf, occ_model) results = tm_dipole( kin, ne, eri, jac_output, property_options={ "operator_A": dipole, "operator_B": dipole, "coordinates": coord, "transition_dipole_moment": True, }, printoptions={"nroot": 2}, ) To include single excitations in the Linear Response, you can use the corresponding :py:class:`~pybest.properties.lr_pccd_s.LRpCCDS` class. .. code-block:: python 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, jac_output, property_options={ "operator_A": dipole, "operator_B": dipole, "coordinates": 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 Jacobian matrix. The results are returned as an :py:class:`~pybest.io.iodata.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 :py:class:`~pybest.io.iodata.IOData` container contains the following attributes. :orb_a: Molecular orbitals from the pCCD reference :e_ref: Total reference energy (pCCD or pCCDS) :e_ee: Excitation energies (first entry is 0.0) :civ_ee: Eigenvectors (CI-like) of the LR Jacobian :t_dm: Transition dipole moments (if requested) Keyword Arguments ================= The LR modules support the following keyword arguments: General Options --------------- :nroot: (int) Number of excited states to compute (default: 1) :threshold: (float) Threshold for printing CI vector contributions (default: 0.1) :dump_cache: (bool) Cache intermediate matrices on disk (default: True for large systems) Davidson Options ---------------- :tolerance: (float) Convergence threshold for excitation energies (default: 1e-6) :tolerancev: (float) Threshold for CI eigenvectors (default: 1e-4) :maxiter: (int) Maximum number of iterations (default: 200) :nguessv: (int) Number of initial guess vectors (default: ``(nroot - 1)*4 + 1``) :maxvectors: (int) Maximum number of vectors in Davidson space (default: ``(nroot - 1)*10``) :todisk: (bool) Write intermediates to disk (default: False) .. _user_properties_lr_examples: Example Usage ============= The example files for obtaining transition dipole moments using the linear response method are collected in the directory ``data/examples/properties``. 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 JacobianpCCD calculation, and an LR-pCCD calculation with a pCCD reference function on the water molecule using the cc-pVDZ basis set. .. literalinclude:: ../src/pybest/data/examples/properties/water_lr_tdm_pccd.py :caption: data/examples/properties/water_lr_tdm_pccd.py :lines: 3- Example for LR-pCCD+S 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 JacobianpCCD+S calculation, and an LR-pCCD+S calculation with a pCCD reference function on the water molecule using the cc-pVDZ basis set. .. literalinclude:: ../src/pybest/data/examples/properties/water_lr_tdm_pccd_s.py :caption: data/examples/properties/water_lr_tdm_pccd_s.py :lines: 3-