..
    : 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 <http://www.gnu.org/licenses/>
    :
    : --

.. _dipole:

Computing the electric dipole moment
####################################

PyBEST provides a method-independent wrapper to calculate the electric dipole
moment of some electronic wave function. The code snippet below shows how to
determine the electric dipole moment for some electronic structure method whose
results are stored in the :py:class:`~pybest.io.iodata.IOData` container
``method_output``. The dipole moment integrals are stored in the ``dipole`` instance
as a ``tuple`` (see :ref:`user_molecularham_dipole` for more details).

.. code-block:: python

    dipole_moment = compute_dipole_moment(dipole, method_output)

The above function returns each component of the dipole moment in the order
{x,y,z}.

.. note::

    In order to calculate the dipole moment, the 1-RDM has to be available and
    stored in the :py:class:`~pybest.io.iodata.IOData` container that is passed
    to the method.

In general, the :py:func:`~pybest.wrapper.multipole.compute_dipole` function
allows you to calculate the dipole moment for a 1-RDM expressed in the atomic
**or** molecular orbital basis. This functionality can be steered using the
``mo`` keyword argument. In total, 3 keyword arguments can be set

    :mo:        (boolean) if True, the 1-RDM is transformed back to the atomic orbital
                 basis (default False)
    :scale:     (float) the scaling factor for the 1-RDM (default 2.0). All
                 internal 1-RDMs in PyBEST are stored for one set of orbitals.
                 Thus, even in the restricted case, only one spin-component is
                 saved. The total 1-RDM (spin-free) is obtained by multiplying
                 by a factor of 2.0
    :name:      (str) if provided, the ``name`` is printed in the header of the
                 std output.


1-RDMs expressed in the atomic orbital basis
============================================

PyBEST stores the following 1-RDMs in the atomic orbital basis:

    - The RHF/UHF 1-RDM

The code snippet below shows how to determine the dipole moment from an RHF
wave function,

.. literalinclude:: ../src/pybest/data/examples/multipole/dipole_water_scf.py
    :lines: 51-58


1-RDMs expressed in the molecular orbital basis
===============================================

All 1-RDMs determined in post-Hartree-Fock methods are stored in the molecular
orbital basis. Thus, they need to be transformed back to the atomic orbital
basis prior to the calculation of the multipole integrals. This back-transformation
is performed automatically if the **keyword argument** ``mo`` is set to ``True``.

The current version of PyBEST supports 1-RDMs for

    - MP2 (relaxed and unrelaxed)
    - OOpCCD
    - LCC (LCCD, LCCSD, pCCD-LCCD, pCCD-LCCSD)

The code snippet below shows how to determine the dipole moment from an OOpCCD
wave function assuming the OOpCCD results are stored in ``oopccd_output``,

.. literalinclude:: ../src/pybest/data/examples/multipole/dipole_water_pccd.py
    :lines: 70-73

The main difference is that we have to set the keyword argument ``mo`` to ``True``.


Example Python scripts
======================

The dipole moment of the water molecule determined from a restricted Hartree-Fock calculation
---------------------------------------------------------------------------------------------

This is a basic example of how to calculate the dipole moment in PyBEST.
This script performs a restricted Hartree-Fock calculation for the water molecule using
the cc-pVDZ basis set.

.. literalinclude:: ../src/pybest/data/examples/multipole/dipole_water_scf.py
    :caption: data/examples/multipole/dipole_water_scf.py
    :lines: 3-


The dipole moment of the water molecule determined from a restricted pCCD calculation
-------------------------------------------------------------------------------------

This is a basic example of how to calculate the dipole moment in PyBEST with the RpCCD method.
This script performs a restricted Hartree-Fock calculation followed by an
orbital-optimized pCCD calculation for the water molecule using
the cc-pVDZ basis set.

.. literalinclude:: ../src/pybest/data/examples/multipole/dipole_water_pccd.py
    :caption: data/examples/multipole/dipole_water_pccd.py
    :lines: 3-73

The dipole moment of the water molecule from RpCCD-LCCSD calculation
--------------------------------------------------------------------

The example code given below shows how to determine the dipole moment from a pCCD-LCCSD wave function.
The script details performing a restricted Hatree-Fock calculation followed by an orbital-optimized
pCCD calculation and then an LCCSD calculation on top of that. Notice that :math:`\Lambda`-amplitudes
are optimized by setting the keyword argument ``l=True`` in the LCC function call.

.. literalinclude:: ../src/pybest/data/examples/multipole/dipole_water_pccd_lccsd.py
    :caption: data/examples/multipole/dipole_water_pccd_lccsd.py
    :lines: 3-