04_CorrectnessCheck.Rmd

---
jupyter:
  jupytext:
    formats: ipynb,Rmd
    text_representation:
      extension: .Rmd
      format_name: rmarkdown
      format_version: '1.2'
      jupytext_version: 1.9.1
  kernelspec:
    display_name: Python 3
    language: python
    name: python3
---

```{python}
# Visualization
import matplotlib.pyplot as plt
from tqdm import tqdm
from time import time
from copy import deepcopy
import matplotlib as matplotlib
from MPS_QFT.helper import get_figsize

# Numpy
import numpy as np
from numpy import linalg as LA

# importing Qiskit
from qiskit import QuantumCircuit, execute, Aer
from qiskit.extensions import Initialize
from qiskit.circuit import library as lb
from qiskit.quantum_info import random_statevector
# %config InlineBackend.figure_format = 'svg' # Makes the images look nice

# Tensor networks
import quimb as quimb
from ncon import ncon

from MPS_QFT.helper import print_state, right_contract, left_contract, reverse_qiskit, qiskit_get_statevect
from MPS_QFT.helper import to_full_MPS, to_dense, to_approx_MPS, Fidelity

from MPS_QFT.manual import apply_two_qubit_gate_full, max_bond_dimension, apply_two_qubit_gate, apply_one_qubit_gate

from MPS_QFT.gates import CPHASE, cphase_swap_qiskit, cphase_swap_quimb

from MPS_QFT.circuit import qft_circuit_qiskit, qft_circuit_swap_full, qft_circuit_swap_approx, circ_data, MPS_circ
from MPS_QFT.circuit import GHZ_qiskit
```

```{python}
# Plotting LateX figures
LateX = False
if LateX:
    matplotlib.use("pgf")
    matplotlib.rcParams.update({ 
        "pgf.texsystem": "pdflatex",
        'font.family': 'serif',
        'text.usetex': True,
        'pgf.rcfonts': False,
        "pgf.preamble": [ r"\usepackage[utf8]{inputenc}" ]
    })

    #Font size configuration
    SMALL_SIZE = 8
    MEDIUM_SIZE = 10
    BIGGER_SIZE = 11
    BIGGEST_SIZE = 12

    #All sizes are customizable here
    plt.rc('font', size=SMALL_SIZE)          # controls default text sizes
    plt.rc('axes', titlesize=BIGGER_SIZE)   # fontsize of the axes title
    plt.rc('axes', labelsize=MEDIUM_SIZE)    # fontsize of the x and y labels
    plt.rc('xtick', labelsize=SMALL_SIZE)    # fontsize of the tick labels
    plt.rc('ytick', labelsize=SMALL_SIZE)    # fontsize of the tick labels
    plt.rc('legend', fontsize=MEDIUM_SIZE)    # legend fontsize
    plt.rc('figure', titlesize=BIGGEST_SIZE)  # fontsize of the figure title
    #plt.rcParams['axes.facecolor'] = 'white'
```

## Tests for approximation error
We take two states, one "simple" with low entanglement (e.g. GHZ), and one random (high entanglement). We then convert them to MPS with a given bond dimension $\chi$, and plot the approximation fidelity as a function of the max bond dimension.
We will use the `_q` notation to indicate the quimb version and `_m` for the manual one.

We then apply a quantum gate (e.g. CNOT) to both states, and do the same plot. Apply the QFT to both states, and repeat again the same plot. 

```{python}
class TestingData():
    """
        Class to eliminate global variables, keeping tidy the informations constant along the testing
    """
    
    def __init__(self, N, d=2):
        """
        Parameters
        ----------
        N: int
            Number of degrees of freedom
        d: int
            Local dimension of dof. If not provided is set to 2, i.e. qubits
        """
        # Number of degrees of freedom
        self.N = N
        # Local dimension of dof
        self.d = d
        # Dimension of the Hilbert space
        self.dim = d**N
```

```{python}
# --- Parameters ---
d = 2 # Qubits
N = 5 # Number of qubits
tx = TestingData(N, d=d)
Chis = np.arange(1, tx.d**(tx.N//2)+4)

# --- Trackers ---
fid_GHZ_q = np.zeros((len(Chis), 3))
fid_ran_q = np.zeros((len(Chis), 3))
fid_GHZ_m = np.zeros((len(Chis), 3))
fid_ran_m = np.zeros((len(Chis), 3))

# --- GHZ state creation ---
qc = QuantumCircuit(tx.N)
GHZ_qiskit(qc) # Works in place
GHZ = qiskit_get_statevect(qc)

# --- Random state creation ---
random_state = random_statevector(tx.dim).data 
qc_ran = Initialize(random_state).definition
random_state4MPS = reverse_qiskit(random_state, tx.N) # Put the state in the correct order
```

```{python}
# --- For the initial state --- 
for i, chi in tqdm(enumerate(Chis)):
    # --- MPS transformation for quimb ---
    GHZ_MPS_q = MPS_circ(qc, chi=chi)
    GHZ_MPS_q =  GHZ_MPS_q.to_dense()[:, 0]
    
    random_MPS_q =  quimb.tensor.tensor_1d.MatrixProductState.from_dense( np.matrix(random_state4MPS), 
                                                                         dims=[d for _ in range(tx.N)], max_bond=chi)
    random_MPS_q = random_MPS_q.to_dense()[:, 0]
    
    # --- Fidelity computations for quimb ---
    fid_GHZ_q[i, 0] = Fidelity(GHZ_MPS_q, GHZ)
    fid_ran_q[i, 0] = Fidelity(random_MPS_q, random_state4MPS ) 
    
    # --- MPS transformation for manual ---
    GHZ_MPS_m    = to_approx_MPS(GHZ, tx.N, chi=chi)
    GHZ_MPS_m    = to_dense(GHZ_MPS_m).reshape(1, tx.dim)
    random_MPS_m = to_approx_MPS(random_state4MPS, N, chi=chi)
    random_MPS_m = to_dense(random_MPS_m).reshape(1, tx.dim)
    
    # --- Fidelity computations for manual ---
    fid_GHZ_m[i, 0] = Fidelity(GHZ_MPS_m, GHZ)
    fid_ran_m[i, 0] = Fidelity(random_MPS_m, random_state4MPS) 
```

```{python}
# --- After the applycation of a CNOT ---
qc = QuantumCircuit(tx.N)
GHZ_qiskit(qc) # Works in place
qc.cx(0, 1)
GHZ_1 = qiskit_get_statevect(qc)

qc_ran = Initialize(random_state).definition
qc_ran.cx( 0, 1)
random_state_1 = qiskit_get_statevect(qc_ran)

for i, chi in tqdm(enumerate(Chis)):
    # --- MPS transformation for quimb ---
    GHZ_MPS_q = MPS_circ(qc, chi=chi)
    GHZ_MPS_q =  GHZ_MPS_q.to_dense()[:, 0]
    random_MPS_q =  quimb.tensor.tensor_1d.MatrixProductState.from_dense( np.matrix(random_state4MPS), 
                                                                         dims=[d for _ in range(tx.N)], max_bond=chi)
    random_MPS_q.gate_(quimb.CNOT(), (0, 1), tags='cx', max_bond=chi, contract='swap+split')
    random_MPS_q = random_MPS_q.to_dense()[:, 0]
    
    # --- Fidelity computations for quimb ---
    fid_GHZ_q[i, 1] = Fidelity(GHZ_MPS_q, GHZ_1)
    fid_ran_q[i, 1] = Fidelity(random_MPS_q, random_state_1) 
    
    # --- MPS transformation for manual ---
    GHZ_MPS_m    = to_approx_MPS(GHZ, tx.N, chi=chi)
    GHZ_MPS_m    = apply_two_qubit_gate( np.array(quimb.CNOT()), 0, GHZ_MPS_m, chi=chi)
    GHZ_MPS_m    = to_dense(GHZ_MPS_m).reshape(1, tx.dim)
    
    random_MPS_m = to_approx_MPS(random_state4MPS, N, chi=chi)
    random_MPS_m = apply_two_qubit_gate( np.array(quimb.CNOT()), 0, random_MPS_m, chi=chi)
    random_MPS_m = to_dense(random_MPS_m).reshape(1, tx.dim)
    
    # --- Fidelity computations for manual ---
    fid_GHZ_m[i, 1] = Fidelity(GHZ_MPS_m, GHZ_1)
    fid_ran_m[i, 1] = Fidelity(random_MPS_m, random_state_1) 
```

```{python}
# --- After the applycation of a QFT ---
qc = QuantumCircuit(tx.N)
GHZ_qiskit(qc) # Works in place
qft_circuit_qiskit(qc, tx.N)
GHZ_2 = qiskit_get_statevect(qc)

qc_ran = Initialize(random_state).definition
qft_circuit_qiskit(qc_ran, tx.N)
random_state_2 = qiskit_get_statevect(qc_ran)

for i, chi in tqdm(enumerate(Chis)):
    # --- MPS transformation for quimb ---
    init = quimb.tensor.tensor_gen.MPS_ghz_state(tx.N)
    qft_circ = QuantumCircuit(tx.N)
    qft_circuit_qiskit(qft_circ, tx.N)
    GHZ_MPS_q = MPS_circ(qft_circ, chi=chi, init_state=init )
    GHZ_MPS_q =  GHZ_MPS_q.to_dense()[:, 0]
    
    random_MPS_q =  quimb.tensor.tensor_1d.MatrixProductState.from_dense( np.matrix(random_state4MPS), 
                                                                         dims=[d for _ in range(tx.N)], max_bond=chi)
    random_MPS_q = MPS_circ(qft_circ, chi=chi, init_state=random_MPS_q)
    random_MPS_q = random_MPS_q.to_dense()[:,0]
    
    # --- Fidelity computations for quimb ---
    fid_GHZ_q[i, 2] = Fidelity(GHZ_MPS_q, GHZ_2)
    fid_ran_q[i, 2] = Fidelity(random_MPS_q, random_state_2) 
    
    # --- MPS transformation for manual ---
    GHZ_MPS_m    = to_approx_MPS(GHZ, tx.N, chi=chi)
    GHZ_MPS_m    = qft_circuit_swap_approx(GHZ_MPS_m, tx.N, chi=chi)
    GHZ_MPS_m    = to_dense(GHZ_MPS_m).reshape(1, tx.dim)
    
    random_MPS_m = to_approx_MPS(random_state4MPS, tx.N, chi=chi)
    random_MPS_m = qft_circuit_swap_approx(random_MPS_m, tx.N, chi=chi)
    random_MPS_m = to_dense(random_MPS_m).reshape(1, tx.dim)
    
    # --- Fidelity computations for manual ---
    fid_GHZ_m[i, 2] = Fidelity(GHZ_MPS_m, GHZ_2)
    fid_ran_m[i, 2] = Fidelity(random_MPS_m, random_state_2) 
```

```{python}
fig, ax = plt.subplots(2, 2, figsize=get_figsize(.95) )
ax = ax.flatten()
colors = plt.cm.rainbow(np.linspace(0,1,3) )

# --- Plotting GHZ in quimb ---
ax[0].plot(Chis, fid_GHZ_q[:, 0], label='Initial GHZ', color=colors[0], linewidth=2)
ax[0].plot(Chis, fid_GHZ_q[:, 1], '--', label='GHZ + CNOT', color=colors[1], linewidth=2)
ax[0].plot(Chis, fid_GHZ_q[:, 2], ':', label='GHZ + QFT', color=colors[2])
ax[0].set_ylabel('Fidelity $|<\\psi|\\phi>|^2$')
ax[0].set_xlabel('Bond dimension  $\\chi$')
ax[0].set_title(f'Quimb, N={tx.N} qubits, GHZ')
ax[0].legend()

# --- Plotting GHZ in manual ---
ax[1].plot(Chis, fid_GHZ_m[:, 0], label='Initial GHZ', color=colors[0], linewidth=2)
ax[1].plot(Chis, fid_GHZ_m[:, 1], '--', label='GHZ + CNOT', color=colors[1], linewidth=2)
ax[1].plot(Chis, fid_GHZ_m[:, 2], ':', label='GHZ + QFT', color=colors[2])
#ax[1].set_ylabel('Fidelity $|<\\psi|\\phi>|^2$')
ax[1].set_xlabel('Bond dimension  $\\chi$')
ax[1].set_title(f'Manual, N={tx.N} qubits, GHZ')
ax[1].legend()

# --- Plotting random in quimb ---
ax[2].plot(Chis, fid_ran_q[:, 0], label='Initial random', color=colors[0], linewidth=2)
ax[2].plot(Chis, fid_ran_q[:, 1], '--', label='Random + CNOT', color=colors[1], linewidth=2)
ax[2].plot(Chis, fid_ran_q[:, 2], ':', label='Random + QFT', color=colors[2])
ax[2].set_ylabel('Fidelity $|<\\psi|\\phi>|^2$')
ax[2].set_xlabel('Bond dimension  $\\chi$')
ax[2].set_title(f'Quimb, N={tx.N} qubits, random')
ax[2].legend()

# --- Plotting random in manual ---
ax[3].plot(Chis, fid_ran_m[:, 0], label='Initial random', color=colors[0], linewidth=2)
ax[3].plot(Chis, fid_ran_m[:, 1], '--', label='Random + CNOT', color=colors[1], linewidth=2)
ax[3].plot(Chis, fid_ran_m[:, 2], ':', label='Random + QFT', color=colors[2])
#ax[3].set_ylabel('Fidelity $|<\\psi|\\phi>|^2$')
ax[3].set_xlabel('Bond dimension  $\\chi$')
ax[3].set_title(f'Manual, N={tx.N} qubits, random')
ax[3].legend()

fig.tight_layout()
if LateX: plt.savefig("Plots/PerfVsChi.pgf")
plt.show()
```

```{python}

```