.. : 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_ipcc_intro: Quick Guide ########### Before reading this section, please read the general introduction mentioned in :ref:`user_posthf_intro`. This part of the Documentation builds upon it. The current version of PyBEST offers ionization potential (IP) calculations with a RpCCD (see :ref:`user_pccd`) reference function using the Equation of Motion (EOM) formalism. The IP module is explained in greater detail below. Supported features ================== The IP module supports spin-restricted orbitals and the ``DenseLinalgFactory`` and ``CholeskyLinalgFactory``. .. _getstartedip: How to: RIP =========== This is a short introduction to the RIP 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 .. note:: In all IP modules, only the spin projection is defined, that is, only the number of excess :math:`\alpha` electrons over :math:`\beta` electrons is defined. .. _qg_ippccd: RHF reference function ---------------------- This version of PyBEST supports single IP-CCD, IP-CCSD, IP-LCCD, and IP-LCCSD with two different excitation operators - IP-CCD (:py:class:`~pybest.ip_eom.xip_rcc.RIPCCD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-CCSD (:py:class:`~pybest.ip_eom.xip_rcc.RIPCCSD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-LCCD (:py:class:`~pybest.ip_eom.xip_rcc.RIPLCCD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-LCCSD (:py:class:`~pybest.ip_eom.xip_rcc.RIPLCCSD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators Complete examples can be found in the following subsections (:ref:`ip_examples`). RpCCD reference function ------------------------ This version of PyBEST supports single IP-pCCD, IP-fpCCD, IP-fpCCSD, IP-fpLCCD, and IP-fpLCCSD and double DIP-pCCD calculations with two different excitation operators - IP-pCCD (:py:class:`~pybest.ip_eom.xip_pccd.RIPpCCD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-fpCCD (:py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-fpCCSD (:py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCSD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-fpLCCD (:py:class:`~pybest.ip_eom.xip_fpcc.RIPfpLCCD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - IP-fpLCCSD (:py:class:`~pybest.ip_eom.xip_pccd.RIPfpLCCSD`): 1. 1 hole (1h) operators 2. 2 hole, 1 particle (2h1p) operators - DIP-pCCD (:py:class:`~pybest.ip_eom.xip_pccd.RDIPpCCD`): 1. 2 hole (2h) operators 2. 3 hole, 1 particle (3h1p) operators Complete examples can be found in the following subsections (:ref:`ip_examples`). IP-CCD/IP-CCSD/IP-LCCD/IP-LCCSD: Doublet and quartet states ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ If you use this module, please cite [galynska2024]_. We assume that you have performed a restricted CC calculation, whose results are stored in the :py:class:`~pybest.io.iodata.IOData` container ``cc_output`` (see :ref:`user_estruct_lcc` for LCC flavors and :ref:`user_estruct_rcc` for conventional CC). The code snippet below shows how to perform a :py:class:`~pybest.ip_eom.xip_rcc.RIPCCSD` calculation. The current version supports the optimization of doublet and quartet states. Both states can be optimized by setting the keyword argument ``alpha`` to ``1`` when creating an instance of :py:class:`~pybest.ip_eom.xip_rcc.RIPCCSD`. The number of targeted states is passed through the ``nroot`` keyword argument in the function call. .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) # Both doublet and quartet states can be targeted ip = RIPCCSD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, cc_output, nroot=3) The results are returned as a :py:class:`~pybest.io.iodata.IOData` container, while all results are saved to the ``pybest-results/checkpoint_IP-EOM-CCSD.h5`` file. Specifically, the :py:class:`~pybest.io.iodata.IOData` container contains (amongst others) the following attributes :e_ip: The ionization energies in :math:`E_h` :civ_ip_alpha: The CI vectors (that is, the eigenvectors) for each state (column) for a given ``alpha`` value The eigenvalues and eigenvectors are stored as numpy arrays, not as instances of the :py:class:`~pybest.linalg.base.LinalgFactory`. Since for one spin projection :math:`S_z`, several spin multiplets are accessible, PyBEST suggests the multiplicity for each targeted state in the output file. By default, :py:class:`~pybest.ip_eom.xip_rcc.RIPCCSD` includes up to 2h1p terms. For ``alpha=1``, the 2h1p terms can be neglected and only the 1h terms are considered during the diagonalization. The number of hole/particle operators is specified by the keyword argument ``nhole``. The following code snippet shows how to perform a :py:class:`~pybest.ip_eom.xip_rcc.RIPCCSD` calculation for 1 unpaired electron including only 1h terms, .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) and 1h terms # Both doublet and quartet states can be targeted ip = RIPCCSD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, cc_output, nroot=3, nhole=1) .. note:: All remaining flavors (RIPCCD, RIPLCCD, and RIPLCCSD) are accessible in a similar way. To choose a different IP-CC flavor, simply change the :py:class:`~pybest.ip_eom.xip_rcc.RIPCCSD` class to the desired CC model, - :py:class:`~pybest.ip_eom.xip_rcc.RIPCCD` for IP-CCD - :py:class:`~pybest.ip_eom.xip_rcc.RIPLCCD` for IP-LCCD - :py:class:`~pybest.ip_eom.xip_rcc.RIPLCCSD` for IP-LCCSD .. note:: A spin-free implementation is also available. It might help to eliminate convergence difficulties in some cases. If set to `True`, only doublet state can be targeted. .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) # Only doublet states can be targeted ip = RIPCCSD(lf, occ_model, alpha=1, spinfree=True) ip_output = ip(kin, ne, eri, cc_output, nroot=3) IP-pCCD: Doublet and quartet states ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ If you use this module, please cite [boguslawski2021]_. We assume that you have performed a restricted pCCD calculation (either with or without orbital optimization), whose results are stored in the :py:class:`~pybest.io.iodata.IOData` container ``pccd_output`` (see :ref:`user_pccd`). The code snippet below shows how to perform a :py:class:`~pybest.ip_eom.xip_pccd.RIPpCCD` calculation. The current version supports the optimization of doublet and quartet states and just quartet states. The former can be optimized by setting the keyword argument ``alpha`` to ``1`` when creating an instance of :py:class:`~pybest.ip_eom.xip_pccd.RIPpCCD`, while the latter are accessible by choosing ``alpha=3``. The number of targeted states is passed through the ``nroot`` keyword argument in the function call. .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) # Both doublet and quartet states can be targeted ip = RIPpCCD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, pccd_output, nroot=3) The results are returned as a :py:class:`~pybest.io.iodata.IOData` container, while all results are saved to the ``pybest-results/checkpoint_IP-EOM-pCCD.h5`` file. Specifically, the :py:class:`~pybest.io.iodata.IOData` container contains (amongst others) the following attributes :e_ip: The ionization energies in :math:`E_h` :civ_ip_alpha: The CI vectors (that is, the eigenvectors) for each state (column) for a given ``alpha`` value The eigenvalues and eigenvectors are stored as numpy arrays, not as instances of the :py:class:`~pybest.linalg.base.LinalgFactory`. Since for one spin projection :math:`S_z`, several spin multiplets are accessible, PyBEST suggests the multiplicity for each targeted state in the output file. By default, :py:class:`~pybest.ip_eom.xip_pccd.RIPpCCD` includes up to 2h1p terms. For ``alpha=1``, the 2h1p terms can be neglected and only the 1h terms are considered during the diagonalization. The number of hole/particle operators is specified by the keyword argument ``nhole``. The following code snippet shows how to perform a :py:class:`~pybest.ip_eom.xip_pccd.RIPpCCD` calculation for 1 unpaired electron including only 1h terms, .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) and 1h terms # Both doublet and quartet states can be targeted ip = RIPpCCD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, pccd_output, nroot=3, nhole=1) .. note:: For three unpaired electrons, only 2h1p terms are supported. Specifying ``nhole=1`` will raise an error. .. note:: A spin-free implementation is also available. It might help to eliminate convergence difficulties in some cases. If set to `True`, only doublet state can be targeted. .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) # Only doublet states can be targeted ip = RIPCCSD(lf, occ_model, alpha=1, spinfree=True) ip_output = ip(kin, ne, eri, cc_output, nroot=3) IP-fpCCD/IP-fpCCSD/IP-fpLCCD/IP-fpLCCSD: Doublet and quartet states ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ If you use this module, please cite [galynska2024]_. We assume that you have performed a restricted frozen-pair (fp)CC calculation, whose results are stored in the :py:class:`~pybest.io.iodata.IOData` container ``cc_output`` (see :ref:`user_estruct_lcc` for LCC flavors on top of pCCD and :ref:`user_estruct_fpcc` for supported fpCC models). The code snippet below shows how to perform a :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCSD` calculation. The current version supports the optimization of doublet and quartet states. Both states can be optimized by setting the keyword argument ``alpha`` to ``1`` when creating an instance of :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCSD`. The number of targeted states is passed through the ``nroot`` keyword argument in the function call. .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) # Both doublet and quartet states can be targeted ip = RIPfpCCSD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, cc_output, nroot=3) The results are returned as a :py:class:`~pybest.io.iodata.IOData` container, while all results are saved to the ``pybest-results/checkpoint_IP-EOM-fpCCSD.h5`` file. Specifically, the :py:class:`~pybest.io.iodata.IOData` container contains (amongst others) the following attributes :e_ip: The ionization energies in :math:`E_h` :civ_ip_alpha: The CI vectors (that is, the eigenvectors) for each state (column) for a given ``alpha`` value The eigenvalues and eigenvectors are stored as numpy arrays, not as instances of the :py:class:`~pybest.linalg.base.LinalgFactory`. Since for one spin projection :math:`S_z`, several spin multiplets are accessible, PyBEST suggests the multiplicity for each targeted state in the output file. By default, :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCSD` includes up to 2h1p terms. For ``alpha=1``, the 2h1p terms can be neglected and only the 1h terms are considered during the diagonalization. The number of hole/particle operators is specified by the keyword argument ``nhole``. The following code snippet shows how to perform a :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCSD` calculation for 1 unpaired electron including only 1h terms, .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) and 1h terms # Both doublet and quartet states can be targeted ip = RIPfpCCSD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, cc_output, nroot=3, nhole=1) .. note:: All remaining flavors (RIPfpCCD, RIPfpLCCD, and RIPfpLCCSD) are accessible in a similar way. To choose a different IP-CC flavor, simply change the :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCSD` class to the desired CC model, - :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpCCD` for IP-fpCCD - :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpLCCD` for IP-fpLCCD - :py:class:`~pybest.ip_eom.xip_fpcc.RIPfpLCCSD` for IP-fpLCCSD .. note:: A spin-free implementation is also available. It might help to eliminate convergence difficulties in some cases. If set to `True`, only doublet state can be targeted. .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.5) # Only doublet states can be targeted ip = RIPfpCCSD(lf, occ_model, alpha=1, spinfree=True) ip_output = ip(kin, ne, eri, cc_output, nroot=3) DIP-pCCD: Singlet, triplet, and quintet states ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ If you use this module, please cite [boguslawski2021]_. We assume that you have performed a restricted pCCD calculation (either with or without orbital optimization), whose results are stored in the :py:class:`~pybest.io.iodata.IOData` container ``pccd_output`` (see :ref:`user_pccd`). The code snippet below shows how to perform a :py:class:`~pybest.ip_eom.xip_pccd.RDIPpCCD` calculation. The current version supports the optimization of singlet, triplet, and quintet states, triplet and quintet states, or just quintet states. The first case can be optimized by setting the keyword argument ``alpha`` to ``0`` when creating an instance of :py:class:`~pybest.ip_eom.xip_pccd.RDIPpCCD`, the second case is accessible by choosing ``alpha=2``, while the last spin multiplicity can be chosen by setting ``alpha=4``. The number of targeted states is passed through the ``nroot`` keyword argument in the function call. .. code-block:: python # Calculate 3 lowest-lying roots for 0 unpaired electron (S_z=0.0) # Singlet, triplet, and quintet states can be targeted ip = RDIPpCCD(lf, occ_model, alpha=0) ip_output = ip(kin, ne, eri, pccd_output, nroot=3) The results are returned as a :py:class:`~pybest.io.iodata.IOData` container, while all results are saved to the ``pybest-results/checkpoint_DIP-EOM-pCCD.h5`` file. Specifically, the :py:class:`~pybest.io.iodata.IOData` container contains (amongst others) the following attributes :e_ip: The ionization energies in :math:`E_h` :civ_ip_alpha: The CI vectors (that is, the eigenvectors) for each state (column) for a given ``alpha`` value The eigenvalues and eigenvectors are stored as numpy arrays, not as instances of the :py:class:`~pybest.linalg.base.LinalgFactory`. Since for one spin projection :math:`S_z`, several spin multiplets are accessible, PyBEST suggests the multiplicity for each targeted state in the output file. By default, :py:class:`~pybest.ip_eom.xip_pccd.RDIPpCCD` includes up to 3h1p terms. For ``alpha=0`` and ``alpha=2``, the 3h1p terms can be neglected so that only the 2h terms are considered during the diagonalization. The number of hole/particle operators is specified by the keyword argument ``nhole``. The following code snippet shows how to perform a :py:class:`~pybest.ip_eom.xip_pccd.RDIPpCCD` calculation for 0 unpaired electron including only 2h terms, .. code-block:: python # Calculate 3 lowest-lying roots for 1 unpaired electron (S_z=0.0) and 2h terms # Singlet and triplet states can be targeted ip = RDIPpCCD(lf, occ_model, alpha=0) ip_output = ip(kin, ne, eri, pccd_output, nroot=3, nhole=2) .. note:: For four unpaired electrons (qunitet states), only 3h1p terms are supported. Specifying ``nhole=2`` will raise an error. .. _ip_fc: Defining a frozen core ---------------------- By default, all core orbitals are frozen. If a frozen core has been selected in the CC reference calculation, the same orbitals have to be frozen in the chosen IP flavor. To freeze some specific (occupied) orbitals, the number of frozen cores has to be specified when an instance of some :py:class:`~pybest.occ_model.OccupationModel` class is created. For instance, if you wish to freeze the first 4 (occupied) orbitals in an IP-pCCD calculation, specify the ``ncore`` argument during the initialization of some occupation model class, .. code-block:: python # Select 4 frozen core orbital #----------------------------- occ_model = AufbauOccModel(gobasis, ncore=4) # Perform IP-CC calculation for 1 unpaired electron (S_z=0.5), # ncore is stored in occ_model #------------------------------------------------------------- ip = RIPpCCD(lf, occ_model, alpha=1) ip_output = ip(kin, ne, eri, pccd_output) This syntax is working for all IP modules mentioned above. .. _ip_restart: Restart options --------------- Restart options are not supported yet.