Emulator
========

Matrix product states
---------------------

.. automodule:: qtealeaves.emulator.mps_simulator
    :members:

Matrix Product State python simulator
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

The Matrix Product States (MPS) are an efficient way of representing a quantum state of 1-dimensional systems. The quantum system
is represented as a tensor network, where each tensor represent a single degree of freedom of dimension **local_dim**. The *entanglement* of the system is
encoded in the links between the tensors, and is bounded by the **maximum bond dimension** :math:`\chi`. The system is evolved by
applying operators to the tensors, which are represented as tensors of dimension *(local_dim x local_dim)*, if they are one-site operators (:py:func:`apply_one_site_operator`), or
as tensors of dimension *(local_dim x local_dim x local_dim x local_dim)*, if they are two-site operators (:py:func:`apply_two_site_operator`).

The entanglement constraint is imposed when we apply the two-site operators. Indeed, to come back to the MPS structure we apply a Singular Value Decomposition (SVD)
and apply a truncation on the singular values. We truncate the singular values :math:`s=\{s_1, s_2, \dots, s_n\}` with :math:`s_1\geq s_2\geq\dots\geq s_n` such that:

.. math::
    \begin{cases}
        s_i \mbox{ is truncated if } & i>\chi \\
        s_i \mbox{ is truncated if } & \frac{s_i}{s_1}\leq \epsilon
    \end{cases}

where :math:`\epsilon` is called **cut ratio** and is usually around :math:`\sim 10^{-9}`.

After the simulation we can perform the following measures:
- Projective measurement of the system, which use a copy of the system and thus allow multiple measurements and keeping the state. Uses :py:func:`meas_projective`.
- Measurement of local one-site operators along the full MPS, using :py:func:`meas_local`.
- Measurements of tensor product operators, which are tensor product of single-site operators. Uses :py:func:`meas_tensor_product`.
- Measurements of weighted sum of tensor product operators, using :py:func:`meas_weighted_sum`.
- The bond  entanglement entropy of the system along each bond, using :py:func:`meas_eantanglement`. This computation adds a minimum overhead, since it is already saved during the simulation and requires only a sum.
For further informations on measuremens refer to :doc:`observables`.

Example
~~~~~~~

Let us prepare an example to understand how to use the simulator. In this example we will use *qubits*, which means we fix the local dimension to 2.
First, we need to decide the number of qubits we are interested in and the convergence parameters. We will only use 2 qubits and build a GHZ state.
Remember that an MPS is always initialized in the **Vacuum** state :math:`|00\dots 0\rangle`, even though you can re-initialize it later on.
(For further informations on the initialization see :py:func:`from_statevector`, :py:func:`from_tensor_list`).

.. code-block::python
    import numpy as np
    from qtealeaves.mps_simulator import MPS
    from tn_py_simulator import TNConvergenceParameters

    Nqubits = 2
    local_dim = 2
    convergence_parameters = TNConvergenceParameters(max_bond_dimension=2, cut_ratio=1e-9) # There are default values, but we fix them
    GHZ_mps = MPS(Nqubits, convergence_parameters, local_dim=local_dim)

Then, we have to evolve the state using one and two-site operators. In particular, we will use the Hadamard and Controlled NOT gates.
We take the transposition of the gates to follow the notation of linear algebra, where the operators are applied to the left of the state,
and so we can multiply the operators :math:`O_i` as follows:

.. math::
    |\psi(t+1)\rangle = O_i |\psi(t)\rangle.

Lets write the code:

.. code-block::python
    # Hadamard
    HH = np.array( [[1, 1], [1, -1]] ).transpose()/np.sqrt(2)
    # Controlled NOT
    CNOT = np.array( [[1 , 0, 0, 0],
                    [0, 1, 0, 0],
                    [0, 0, 0, 1],
                    [0, 0, 1, 0]
            ]).transpose()

    # Apply the Hadamard to qubit 0
    GHZ_mps.apply_one_site_operator(HH, 0)
    # Apply the CNOT operator
    # Notice that the apply_two_site_operator requires a tensor of shape
    # (dim_local, dim_local, dim_local, dim_local), and so we need to reshape
    GHZ_mps.apply_two_site_operator(CNOT.reshape(2, 2, 2, 2), pos=0, swap=False)

It is interesting to notice that we used a parameter called `swap` and set it to `False`. The :py:func:`apply_two_site_operator` method
always applies the gates on position `pos, pos+1`. If the **control** qubit of the gate is on index **pos** then the swap parameter should
be set to `False`, while if the **control** is in **pos+1** it should be set to `True`. This function performs the swapping by permuting the
tensor.
After evolving the state, we can perform measuremnts or obtain the real statevector for debugging:

.. code-block::python
    proj_meas = GHZ_mps.meas_projective(nmeas=1024)
    statevector = GHZ_mps.to_stavector()
    for keys, values in zip( proj_meas.keys(), proj_meas.values() ):
        print(keys, ':', values)

    >>> 00 : 502
    >>> 11 : 522

As expected, we got a measurement compatible with the state :math:`|GHZ\rangle=\frac{1}{\sqrt{2}}(|00\rangle+|11\rangle)`.

Finally, remember that you can access the tensors of the MPS class as you would access the elements of a list, or concatenate
two mps using the :py:method:`kron` method.

.. code-block::python
    # first tensor
    first_tens = GHZ_mps[0]
    # last tensor
    last_tens = GHZ_mps[-1]
    # full list of tensors
    tensors = GHZ_mps[:]
    # concatenation
    concMPS = GHZ_mps.kron( GHZ_mps )

List of overloaded operators
~~~~~~~~~~~~~~~~~~~~~~~~~~~~

We report here a list of the overloaded operators and their actions, for an easy guide:

- `+`, add two MPS states. If the first is in the state :math:`|00\rangle` and the second in :math:`|11\rangle`, then their sum is :math:`|GHZ\rangle=\frac{1}{\sqrt{2}}(|00\rangle+|11\rangle)`
- `+=`, as before but acts in place
- `*`, multiply the gauge/orthogonality center by a scalar number
- `*=`, as before but acts in place
- `/`, divide the gauge/orthogonality center by a scalar number
- `/=`, as before but acts in place
- `@`, contract two mps along all the open indexes. The left part is complex conjugate, which means `a @ b` is equivalent to :math:`\langle a | b \rangle`


Tree tensor networks
--------------------

.. autoclass:: qtealeaves.emulator.TTN
   :members:


Locally purified tensor networks
--------------------------------

.. autoclass:: qtealeaves.emulator.LPTN
   :members:


Tree tensor operators
---------------------

.. autoclass:: qtealeaves.emulator.TTO
   :members:


Augmented tree tensor networks
------------------------------

.. autoclass:: qtealeaves.emulator.ATTN
   :members:


MPI version of matrix product states
------------------------------------

.. autoclass:: qtealeaves.emulator.MPIMPS
   :members:


State vectors
-------------

.. autoclass:: qtealeaves.emulator.StateVector
   :members:


Tensor network nodes
--------------------

.. autoclass:: qtealeaves.emulator.TNnode
   :members:


Abstract tensor networks
------------------------

.. autoclass:: qtealeaves.abstracttns.abstract_tn._AbstractTN
   :members:

.. autoclass:: qtealeaves.abstracttns._AbstractMatrixTN
   :members:

.. autoclass:: qtealeaves.abstracttns.postprocess_statedict
   :members:


Setup unitaries before projective measurements
----------------------------------------------

.. autoclass:: qtealeaves.emulator.UnitarySetupProjMeas
   :members:
