.. : 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_rsfcc_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 Reversed Spin Flip (RSF) calculations with RCCD and RCCSD (see :ref:`user_estruct_rcc`), RLCCD and RLCCSD (see :ref:`user_estruct_lcc`), RfpCCD and RfpCCSD (see :ref:`user_estruct_fpcc`), and RSFfpLCCD and RSFfpLCCSD (see :ref:`user_estruct_fpcc`) reference functions using the Equation of Motion (EOM) formalism. The RSF module is explained in greater detail below. Supported features ================== The RSF module supports spin-restricted orbitals and the ``DenseLinalgFactory`` and ``CholeskyLinalgFactory``. .. _getstartedrsf: How to: RSF =========== This is a short introduction to the RSF 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 RSF modules, only the spin projection is defined, that is, only the number of excess :math:`\alpha` electrons over :math:`\beta` electrons is defined. .. _qg_rsfpccd: RCCD-type reference function ---------------------------- This version of PyBEST supports RSF-CCD, RSF-LCCD, RSF-RfpCCD and RSF-RfpLCCD calculations for two different states - RSF-CCD (:py:class:`~pybest.rsf_eom.xrsf_cc.RSFCCD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``To Be Announced (TBA)`` - RSF-LCCD (:py:class:`~pybest.rsf_eom.xrsf_cc.RSFLCCD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` - RSF-RfpCCD (:py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpCCD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` - RSF-RfpLCCD (:py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpLCCD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` Complete examples can be found in the following subsections (:ref:`rsf_examples`). RCCSD-type reference function ----------------------------- This version of PyBEST supports RSF-CCSD, RSF-LCCSD, RSF-RfpCCSD and RSF-RfpLCCSD calculations for two different states - RSF-CCSD (:py:class:`~pybest.rsf_eom.xrsf_cc.RSFCCSD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` - RSF-LCCSD (:py:class:`~pybest.rsf_eom.xrsf_cc.RSFLCCSD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` - RSF-RfpCCSD (:py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpCCSD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` - RSF-RfpLCCSD (:py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpLCCSD`): 1. quintet state (`S_z=2`) with `Ms=2` 2. triplet state (`S_z=1`) with `Ms=1` ``(TBA)`` Complete examples can be found in the following subsections (:ref:`rsf_examples`). RSF-CCD/RSF-CCSD/RSF-LCCD/RSF-LCCSD Triplet and Quintet states ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ If you use this module, please cite [sujkowski2024]_. 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.rsf_eom.xrsf_cc.RSFCCSD` calculation. The current version supports the optimization of triplet (TBA) and quintet states. The state can be optimized by setting the keyword argument ``alpha``: ``2`` for triplet and ``4`` for quintet when creating an instance of :py:class:`~pybest.rsf_eom.xrsf_cc.RSFCCSD`. The number of targeted states is passed through the ``nroot`` keyword argument in the function call. .. code-block:: python # Do RSF-CCSD calculation for 4 unpaired electrons (quintet) rsf = RSFCCSD(lf, occ_model, alpha=4) rsf_output = rsf(kin, ne, eri, ccsd_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_RSF-EOM-CCSD.h5`` file. Specifically, the :py:class:`~pybest.io.iodata.IOData` container contains (amongst others) the following attributes :e_ee: The excitation energies (including the ground state, which equals 0.0) :civ_ee: The CI vector (eigenvectors) for each state (column) with the ground state as the first column vector 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. .. note:: All remaining flavors (RSFCCD, RSFLCCD, and RSFLCCSD) are accessible in a similar way. To choose a different RSF-CC flavor, simply change the :py:class:`~pybest.rsf_eom.xrsf_cc.RSFCCSD` class to the desired CC model, - :py:class:`~pybest.rsf_eom.xrsf_cc.RSFCCD` for RSF-CCD - :py:class:`~pybest.rsf_eom.xrsf_cc.RSFLCCD` for RSF-LCCD - :py:class:`~pybest.rsf_eom.xrsf_cc.RSFLCCSD` for RSF-LCCSD RSF-fpCCD/RSF-fpCCSD/RSF-fpLCCD/RSF-fpLCCSD Triplet and Quintet states ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ If you use this module, please cite [sujkowski2024]_. We assume that you have performed a restricted fpCC 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.rsf_eom.xrsf_cc.RSFfpCCSD` calculation. The current version supports the optimization of triplet (TBA) and quintet states. The state can be optimized by setting the keyword argument ``alpha``: ``2`` for triplet and ``4`` for quintet when creating an instance of :py:class:`~pybest.rsf_eom.xrsf_cc.RSFfpCCSD`. The number of targeted states is passed through the ``nroot`` keyword argument in the function call. .. code-block:: python # Do RSF-fpCCSD calculation for 4 unpaired electrons (quintet) rsf_fp = RSFfpCCSD(lf, occ_model, alpha=4) rsf_fp_output = rsf_fp(kin, ne, eri, fpccsd_output, nroot=4) The results are returned as a :py:class:`~pybest.io.iodata.IOData` container, while all results are saved to the ``pybest-results/checkpoint_RSF-EOM-fpCCSD.h5`` file. Specifically, the :py:class:`~pybest.io.iodata.IOData` container contains (amongst others) the following attributes :e_ee: The excitation energies (including the ground state, which equals 0.0) :civ_ee: The CI vector (eigenvectors) for each state (column) with the ground state as the first column vector 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. .. note:: All remaining flavors (RSFfpCCD, RSFfpLCCD, and RSFfpLCCSD) are accessible in a similar way. To choose a different RSF-fpCC flavor, simply change the :py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpCCSD` class to the desired CC model, - :py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpCCD` for RSF-fpCCD - :py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpLCCD` for RSF-fpLCCD - :py:class:`~pybest.rsf_eom.xrsf_fpcc.RSFfpLCCSD` for RSF-fpLCCSD