Manybody calculations
=====================================

Manybody calculations aim to compute interaction energies and to extrapolate total system energies
and properties by computing the energies/properties on collections of fragments that compose a larger
system.

For more details about manybody calculations in general, see blah.

Manybody computations in QCArchive are implemented as a :ref:`glossary_service`
and computed using the `QCManyBody <https://github.com/MolSSI/QCManybody>`_
package. The individual calculations themselves can be run from any of the QM packages available
in QCArchive.

See also `the QCManybody article <https://doi.org/10.1063/5.0231843>`_

.. note::

     The "initial molecules" for the manybody calculations must contain more than one fragment. It's these
     fragments that are used to construct the clusters for the manybody calculations.


.. _manybody_record:

Manybody Record
-------------------

The :class:`~qcportal.manybody.record_models.ManybodyRecord` object contains all the fields of the
:doc:`base record class <base>`, but also contains manybody-specific fields:

- ``initial_molecule`` - The initial (fragmented molecule) that the manybody calculation was run on
- ``specification`` - The level of theory and other options for running the calculation (see below)
- ``clusters`` - List of :class:`~qcportal.manybody.record_models.ManybodyCluster` objects that represent the
  individual (child) calculations (that is, the calculations run on all the monomers, dimers, etc).

The :class:`~qcportal.manybody.record_models.ManybodyCluster` objects contain:

- ``molecule_id`` and ``molecule`` - the actual molecule this child calculation used. This molecule is a subset
  of the fragments of the input molecule.
- ``fragments`` - Indices of the fragments this molecule/calculation represents
- ``basis`` - The fragment basis used in this molecule. The molecule may contain ghost atoms for
  some centers to deal with basis-set superposition error, and the ghosted fragments
  represent additions to the basis. For example, ``fragments=[4]`` and ``basis=[4,8]`` means
  the calculation is run on fragment 4, but includes basis sets from both fragments 4 and 8.
  Note that this is fragments, not atom indices.
- ``mc_level`` - Label of the model chemistry (ie, 'psi4/b3lyp/6-31g*'). This is automatically generated by QCArchive.
- ``singlepoint_id`` and ``singlepoint_record`` - The :ref:`singlepoint record <singlepoint_record>` of the calculation.

The full calculated results of the manybody calculation itself are stored in the `properties` field of the record.


.. _manybody_specification:

Manybody Specification
-------------------------------------------

The :ref:`glossary_specification` for a manybody record is a :class:`~qcportal.manybody.record_models.ManybodySpecification`.
The fields of that are:

- ``program`` - Program to run the manybody calculation itself. Only value available right now is ``qcmanybody``
- ``levels`` - Dictionary of the different singlpoint specifications for different levels of the expansion (see below)
- ``bsse_correction`` - Correction for the basis set superposition error. This is a list that can contain multiple values,
  enabling calculating with different levels of correction simulataneously. Values are:

  - ``nocp`` for no correction
  - ``cp`` for counterpoise correction
  - ``vmfc`` for Valiron-Mayer Function Counterpoise

- ``keywords`` - Additional keywords for running the calculation. Currently, one keyword is supported:

  - `return_total_data` - Computes/Returns the total energy/gradient, etc, of the system. If false, only the interaction
    energy is returned. Default is false.


- ``protocols`` - Additional keywords controlling the return of information. This is reserved for future use.

The ``levels`` parameter is a dictionary of levels of expansion to a specification describing
how that level should be computed. The keys are integers representing the level (1 = monomers, 2 = dimers, etc)
and the value is a :ref:`singlepoint specification/qc specification <singlepoint_specification>`. A string
"supersystem" key can be used for the supersystem calculation.


.. dropdown:: Compute the interaction energy of dimers using one level of theory

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        from qcportal.manybody import ManybodySpecification

        ManybodySpecification(
            program="qcmanybody",
            bsse_correction=["nocp"],
            levels={
                1: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                ),
                2: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                ),
            },
        )


.. dropdown:: Interaction energy as above, but adding CP correction

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        from qcportal.manybody import ManybodySpecification

        ManybodySpecification(
            program="qcmanybody",
            bsse_correction=["nocp", "cp"],
            levels={
                1: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                ),
                2: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                ),
            },
        )


.. dropdown:: Include total energy for a 4-fragment system

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        from qcportal.manybody import ManybodySpecification

        ManybodySpecification(
            program="qcmanybody",
            bsse_correction=["cp"],
            levels={
                1: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                ),
                2: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                ),
            },
            keywords={'return_total_data': True}
        )


.. dropdown:: Using different levels of theory

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        from qcportal.manybody import ManybodySpecification

        ManybodySpecification(
            program="qcmanybody",
            bsse_correction=["nocp", "cp"],
            levels={
                1: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="ccsd",
                    basis="cc-pvqz",
                ),
                2: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="ccsd",
                    basis="cc-pvtz",
                ),
                3: QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="ccsd",
                    basis="cc-pvdz",
                ),
            }
        )


.. _manybody_submission:

Submitting Records
------------------

Manybody records can be submitted using a client via the :meth:`~qcportal.client.PortalClient.add_manybodys` method.
This method takes the following information:

- ``initial_molecules`` - A single molecule or list of molecules to compute
- ``program`` - Always ``qcmanybody``
- ``levels``, ``bsse_correction``, ``keywords`` - Same as in the :ref:`specification <manybody_specification>` above

See :doc:`../record_submission` for more information about other fields.


.. _manybody_dataset:

Manybody Dataset
--------------------

Manybody :ref:`datasets <glossary_dataset>` are collections of manybody records.
:class:`Entries <qcportal.manybody.dataset_models.ManybodyDatasetEntry>` contain the initial (fragmented)
molecule.

The :class:`dataset specifications <qcportal.manybody.dataset_models.ManybodyDatasetSpecification>`
contain a manybody specification (see :ref:`above <manybody_specification>`)


.. _manybody_client_examples:

Client Examples
---------------

When creating the ``levels`` member of the specification or to be passed into the
:meth:`~qcportal.client.PortalClient.add_manybodys` client function, it is often
easiest to define the specification separately and use the resulting object
as values in the ``levels`` dictionary.

.. dropdown:: Obtain a single record by ID

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        r = client.get_manybodys(123)

.. dropdown:: Obtain multiple records by ID

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        r_lst = client.get_manybodys([123, 456])

.. dropdown:: Obtain multiple records by ID, ignoring missing records

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        r_lst = client.get_manybodys([123, 456, 789], missing_ok=True)

.. dropdown:: Include all data for a record during initial fetch

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        r_lst = client.get_manybodys([123, 456], include=['**'])

.. dropdown:: Query manybody records by program, method, basis

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        # NOTE - we are querying by the program actually used in the singlepoint calculations here
        #        not the overarching manybody program

        r_iter = client.query_manybodys(qc_program='psi4', qc_method='b3lyp', qc_basis='def2-svp')
        for r in r_iter:
            print(r.id)

.. dropdown:: Query manybody records by program and when the record was created, include all data

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        # NOTE - we are querying by the program actually used in the singlepoint calculations here
        #        not the overarching manybody program

        r_iter = client.query_manybodys(program='psi4',
                                        created_after='2024-03-21 12:34:56',
                                        include=['**'])
                                        limit=50)
        for r in r_iter:
            print(r.id)

.. dropdown:: Create some fragmented molecules

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        water_dimer = Molecule(symbols=['O', 'H', 'H', 'O', 'H', 'H'],
                               geometry=[2.81211079,  0.1255717,   0.0,
                                         3.48216662, -1.5543998,   0.0,
                                         1.00578203, -0.1092573,   0.0,
                                        -2.68215279, -0.12325075,  0.0,
                                        -3.27523823,  0.81341093,  1.43347254,
                                        -3.27523823,  0.81341093, -1.43347254],
                               fragments=[[0, 1, 2], [3, 4, 5]],
                              )

        water_stacked = Molecule(symbols=['O', 'H', 'H', 'O', 'H', 'H', 'O', 'H', 'H', 'O', 'H', 'H'],
                                 geometry=[2.81211079,  0.1255717,  0.0,
                                           3.48216662, -1.5543998,  0.0,
                                           1.00578203, -0.1092573,  0.0,
                                           2.81211079,  0.1255717,  4.0,
                                           3.48216662, -1.5543998,  4.0,
                                           1.00578203, -0.1092573,  4.0,
                                           2.81211079,  0.1255717,  8.0,
                                           3.48216662, -1.5543998,  8.0,
                                           1.00578203, -0.1092573,  8.0,
                                           2.81211079,  0.1255717, 12.0,
                                           3.48216662, -1.5543998, 12.0,
                                           1.00578203, -0.1092573, 12.0
                                 ],
                                 fragments=[[0, 1, 2], [3, 4, 5], [6, 7, 8], [9, 10, 11]],
                              )


.. dropdown:: Add a manybody record - compute the interaction energy of a dimer, no BSSE correction

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        qc_spec = QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                  )

        meta, ids = client.add_manybodys([water_dimer],
                                          program='qcmanybody',
                                          levels={
                                            1: qc_spec,
                                            2: qc_spec,
                                          },
                                          bsse_correction=['nocp'],
                                          keywords={}
                                         )


.. dropdown:: Add a manybody record - compute the interaction energy of a dimer, CP-corrected

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        qc_spec = QCSpecification(
                    program="psi4",
                    driver="energy",
                    method="b3lyp",
                    basis="6-31G*",
                  )

        meta, ids = client.add_manybodys([water_dimer],
                                          program='qcmanybody',
                                          levels={
                                            1: qc_spec,
                                            2: qc_spec,
                                          },
                                          bsse_correction=['cp'],
                                          keywords={}
                                         )

.. dropdown:: Add a manybody record - using different levels of theory

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        qc_spec_1 = QCSpecification(
                      program="psi4",
                      driver="energy",
                      method="b3lyp",
                      basis="def2-tzvp",
                    )

        qc_spec_2 = QCSpecification(
                      program="psi4",
                      driver="energy",
                      method="b3lyp",
                      basis="def2-svp",
                    )

        meta, ids = client.add_manybodys([water_dimer],
                                          program='qcmanybody',
                                          levels={
                                            1: qc_spec_1,
                                            2: qc_spec_2,
                                          },
                                          bsse_correction=['nocp'],
                                          keywords={}
                                         )


.. _manybody_dataset_examples:

Dataset Examples
----------------

.. dropdown:: Create a manybody dataset with default options

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        ds = client.add_dataset(
                 "manybody",
                 "Dataset Name",
                 "An example of a manybody dataset"
        )

.. dropdown:: Add a single entry to a manybody dataset

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        water_dimer = Molecule(symbols=['O', 'H', 'H', 'O', 'H', 'H'],
                               geometry=[2.81211079,  0.1255717,   0.0,
                                         3.48216662, -1.5543998,   0.0,
                                         1.00578203, -0.1092573,   0.0,
                                        -2.68215279, -0.12325075,  0.0,
                                        -3.27523823,  0.81341093,  1.43347254,
                                        -3.27523823,  0.81341093, -1.43347254],
                               fragments=[[0, 1, 2], [3, 4, 5]],
                              )

        ds.add_entry("water_dimer", water_dimer)

.. dropdown:: Add many entries to a manybody dataset

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        from qcportal.manybody import ManybodyDatasetEntry

        # Construct a list of entries to add somehow
        ent_1 = ManybodyDatasetEntry(name="water_dimer", initial_molecule=water_dimer)
        ent_2 = ManybodyDatasetEntry(name="water_stacked", initial_molecule=water_stacked)

        new_entries = [ent_1, ent_2]

        # Efficiently add all entries in a single call
        ds.add_entries(new_entries)

.. dropdown:: Add a specification to a manybody dataset

  .. tab-set::

    .. tab-item:: PYTHON

      .. code-block:: py3

        from qcportal.manybody import ManybodySpecification

        spec = ManybodySpecification(
                   program="qcmanybody",
                   bsse_correction=["nocp", "cp"],
                   levels={
                       1: QCSpecification(
                           program="psi4",
                           driver="energy",
                           method="ccsd",
                           basis="cc-pvqz",
                       ),
                       2: QCSpecification(
                           program="psi4",
                           driver="energy",
                           method="ccsd",
                           basis="cc-pvtz",
                       ),
                       3: QCSpecification(
                           program="psi4",
                           driver="energy",
                           method="ccsd",
                           basis="cc-pvdz",
                       ),
                   }
               )

        ds.add_specification("psi4/ccsd/multi", spec)


.. _manybody_api_links:

Manybody QCPortal API
------------------------

* :mod:`Record models <qcportal.manybody.record_models>`
* :mod:`Dataset models <qcportal.manybody.dataset_models>`

* PortalClient methods

  * :meth:`~qcportal.client.PortalClient.add_manybodys`
  * :meth:`~qcportal.client.PortalClient.get_manybodys`
  * :meth:`~qcportal.client.PortalClient.query_manybodys`