Chapter 2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "bf1dd8e0",
   "metadata": {},
   "source": [
    "# The QuTiP framework\n",
    "\n",
    "The QuTiP package can be imported with"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "5a899740",
   "metadata": {},
   "outputs": [],
   "source": [
    "import qutip"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0c0f4d03",
   "metadata": {},
   "source": [
    "It can also be imported with the command `from qutip import *`, that automatically imports all the QuTiP functions. However, here we use the first method, in order to explicitly see the QuTiP functions."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "08dbdf93",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "QuTiP: Quantum Toolbox in Python\n",
      "================================\n",
      "Copyright (c) QuTiP team 2011 and later.\n",
      "Current admin team: Alexander Pitchford, Nathan Shammah, Shahnawaz Ahmed, Neill Lambert, Eric Giguère, Boxi Li, Jake Lishman, Simon Cross and Asier Galicia.\n",
      "Board members: Daniel Burgarth, Robert Johansson, Anton F. Kockum, Franco Nori and Will Zeng.\n",
      "Original developers: R. J. Johansson & P. D. Nation.\n",
      "Previous lead developers: Chris Granade & A. Grimsmo.\n",
      "Currently developed through wide collaboration. See https://github.com/qutip for details.\n",
      "\n",
      "QuTiP Version:      5.0.0\n",
      "Numpy Version:      1.26.4\n",
      "Scipy Version:      1.12.0\n",
      "Cython Version:     3.0.8\n",
      "Matplotlib Version: 3.8.0\n",
      "Python Version:     3.12.1\n",
      "Number of CPUs:     32\n",
      "BLAS Info:          Generic\n",
      "INTEL MKL Ext:      False\n",
      "Platform Info:      Linux (x86_64)\n",
      "Installation path:  /home/alberto/anaconda3/envs/physics/lib/python3.12/site-packages/qutip\n",
      "================================================================================\n",
      "Please cite QuTiP in your publication.\n",
      "================================================================================\n",
      "For your convenience a bibtex reference can be easily generated using `qutip.cite()`\n"
     ]
    }
   ],
   "source": [
    "qutip.about()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3e0cb362",
   "metadata": {},
   "source": [
    "## Quantum Operators\n",
    "\n",
    "Quantum operators play a crucial role in the formulation of quantum mechanics, representing physical observables and operations on quantum states. In QuTiP, operators are represented as `Qobj` instances, just like quantum states. This section introduces the creation and manipulation of quantum operators.\n",
    "\n",
    "### Creating Operators\n",
    "\n",
    "Operators in quantum mechanics can represent measurements, such as position or momentum, and transformations, such as rotation. Let's see how we can define some common operators in QuTiP.\n",
    "\n",
    "#### The Annihilation Operator of the Quantum Harmonic oscillator\n",
    "\n",
    "The harmonic oscillator is a fundamental model in quantum mechanics for understanding various physical systems. Its quantization leads to the concept of creation and annihilation operators, which respectively increase and decrease the energy of the system by one quantum of energy.\n",
    "\n",
    "The annihilation operator, often denoted by $\\hat{a}$, acts on a quantum state to reduce its quantum number. The action of $\\hat{a}$ on a state $|n\\rangle$ is defined as:\n",
    "\n",
    "$$\\hat{a} |n\\rangle = \\sqrt{n} |n-1\\rangle$$\n",
    "\n",
    "Here, $|n\\rangle$ represents a quantum state with $n$ quanta of energy (also known as a Fock state), and $\\sqrt{n}$ is the normalization factor. The matrix representation of the annihilation operator in an $N$-dimensional Hilbert space is given by an upper triangular matrix with the square roots of natural numbers as its off-diagonal elements."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "297c0229",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Annihilation operator (a) for N=7:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [10]], shape=(10, 10), type='oper', dtype=Dia, isherm=False$$\\left(\\begin{array}{cc}0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 1.414 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 1.732 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 2 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 2.236 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 2.449 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 2.646 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 2.828 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 3\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [10]], shape=(10, 10), type='oper', dtype=Dia, isherm=False\n",
       "Qobj data =\n",
       "[[0.         1.         0.         0.         0.         0.\n",
       "  0.         0.         0.         0.        ]\n",
       " [0.         0.         1.41421356 0.         0.         0.\n",
       "  0.         0.         0.         0.        ]\n",
       " [0.         0.         0.         1.73205081 0.         0.\n",
       "  0.         0.         0.         0.        ]\n",
       " [0.         0.         0.         0.         2.         0.\n",
       "  0.         0.         0.         0.        ]\n",
       " [0.         0.         0.         0.         0.         2.23606798\n",
       "  0.         0.         0.         0.        ]\n",
       " [0.         0.         0.         0.         0.         0.\n",
       "  2.44948974 0.         0.         0.        ]\n",
       " [0.         0.         0.         0.         0.         0.\n",
       "  0.         2.64575131 0.         0.        ]\n",
       " [0.         0.         0.         0.         0.         0.\n",
       "  0.         0.         2.82842712 0.        ]\n",
       " [0.         0.         0.         0.         0.         0.\n",
       "  0.         0.         0.         3.        ]\n",
       " [0.         0.         0.         0.         0.         0.\n",
       "  0.         0.         0.         0.        ]]"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define the annihilation operator for N-dimensional Hilbert space\n",
    "N = 10\n",
    "\n",
    "a = qutip.destroy(N)\n",
    "\n",
    "print(\"Annihilation operator (a) for N=7:\")\n",
    "a"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bfde26e3",
   "metadata": {},
   "source": [
    "#### Pauli Matrices\n",
    "\n",
    "The Pauli matrices are fundamental in the study of quantum mechanics, representing the spin operators for a spin-1/2 particle and quantum two-level systems.\n",
    "\n",
    "$$\\sigma_x = \\begin{pmatrix} 0 & 1 \\\\ 1 & 0 \\end{pmatrix}, \\ \\sigma_y = \\begin{pmatrix} 0 & -i \\\\ i & 0 \\end{pmatrix}, \\ \\sigma_z = \\begin{pmatrix} 1 & 0 \\\\ 0 & -1 \\end{pmatrix}$$\n",
    "\n",
    "We can define these matrices in QuTiP as follows:\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "9c4d27f3",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sigma X:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=True$$\\left(\\begin{array}{cc}0 & 1\\\\1 & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=True\n",
       "Qobj data =\n",
       "[[0. 1.]\n",
       " [1. 0.]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "Sigma Y:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=True$$\\left(\\begin{array}{cc}0 & -1j\\\\1j & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=True\n",
       "Qobj data =\n",
       "[[0.+0.j 0.-1.j]\n",
       " [0.+1.j 0.+0.j]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "Sigma Z:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=True$$\\left(\\begin{array}{cc}1 & 0\\\\0 & -1\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=True\n",
       "Qobj data =\n",
       "[[ 1.  0.]\n",
       " [ 0. -1.]]"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sigma_x = qutip.sigmax()\n",
    "sigma_y = qutip.sigmay()\n",
    "sigma_z = qutip.sigmaz()\n",
    "\n",
    "print(\"Sigma X:\")\n",
    "display(sigma_x)\n",
    "print(\"\\n\")\n",
    "print(\"Sigma Y:\")\n",
    "display(sigma_y)\n",
    "print(\"\\n\")\n",
    "print(\"Sigma Z:\")\n",
    "sigma_z"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "037c82c4",
   "metadata": {},
   "source": [
    "### Operator Functions and Operations\n",
    "\n",
    "QuTiP supports various operations on operators, including addition, multiplication (both scalar and matrix), and the commutator. These operations are essential for constructing Hamiltonians, calculating observables, and more.\n",
    "\n",
    "#### Example: Commutator of Pauli Matrices\n",
    "\n",
    "The commutator of two operators $A$ and $B$ is defined as $[A, B] = AB - BA$. Let's calculate the commutator of $\\sigma_x$ and $\\sigma_y$.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "b600c559",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Commutator of Sigma X and Sigma Y:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=False$$\\left(\\begin{array}{cc}2j & 0\\\\0 & -2j\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=CSR, isherm=False\n",
       "Qobj data =\n",
       "[[0.+2.j 0.+0.j]\n",
       " [0.+0.j 0.-2.j]]"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "commutator_xy = qutip.commutator(sigma_x, sigma_y)\n",
    "print(\"Commutator of Sigma X and Sigma Y:\")\n",
    "commutator_xy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "780d6d02",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "commutator_xy == 2j * sigma_z"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9924cf37",
   "metadata": {},
   "source": [
    "## Quantum States\n",
    "\n",
    "Quantum states describe the state of a quantum system. In QuTiP, states are represented again as `Qobj` instances. This section focuses on the representation and manipulation of quantum states.\n",
    "\n",
    "### Fock States\n",
    "\n",
    "The most basic quantum states are the fock states, often denoted as $|n\\rangle$ (with $n \\in \\mathbb{N}$). Let's see how we can create these in QuTiP.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "9c6e2bdf",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "|0> state:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense$$\\left(\\begin{array}{cc}1\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense\n",
       "Qobj data =\n",
       "[[1.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "|1> state:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense$$\\left(\\begin{array}{cc}0\\\\1\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense\n",
       "Qobj data =\n",
       "[[0.]\n",
       " [1.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]\n",
       " [0.]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Creating fock states\n",
    "fock_0 = qutip.fock(N, 0)  # |0>\n",
    "fock_1 = qutip.fock(N, 1)  # |1>\n",
    "\n",
    "print(\"|0> state:\")\n",
    "display(fock_0)\n",
    "print(\"\\n\")\n",
    "print(\"|1> state:\")\n",
    "display(fock_1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cc8b7add",
   "metadata": {},
   "source": [
    "### Superposition States\n",
    "\n",
    "Quantum mechanics allows particles to be in a superposition of states. Let's create a superposition state.\n",
    "\n",
    "$$\n",
    "\\vert \\psi \\rangle = \\frac{1}{\\sqrt{2}} \\left( \\vert 0 \\rangle + \\vert 1 \\rangle \\right)\n",
    "$$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "7db22e3e",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Superposition state:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense$$\\left(\\begin{array}{cc}0.707\\\\0.707\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\\\0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense\n",
       "Qobj data =\n",
       "[[0.70710678]\n",
       " [0.70710678]\n",
       " [0.        ]\n",
       " [0.        ]\n",
       " [0.        ]\n",
       " [0.        ]\n",
       " [0.        ]\n",
       " [0.        ]\n",
       " [0.        ]\n",
       " [0.        ]]"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Creating a superposition state\n",
    "superposition_state = (fock_0 + fock_1).unit()  # Normalize the state\n",
    "\n",
    "print(\"Superposition state:\")\n",
    "superposition_state"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7b554ab1",
   "metadata": {},
   "source": [
    "### Coherent States\n",
    "\n",
    "Coherent states in QuTiP represent quantum states closest to classical waves, defined as\n",
    "$$\n",
    "|\\alpha\\rangle = e^{-|\\alpha|^2/2} \\sum_{n=0}^{\\infty} \\frac{\\alpha^n}{\\sqrt{n!}} |n\\rangle \\, ,\n",
    "$$\n",
    "with minimal uncertainty.\n",
    "\n",
    "The coherent state is an eigenstate of the annihilation operator\n",
    "$$\n",
    "\\hat{a} \\vert \\alpha \\rangle = \\alpha \\vert \\alpha \\rangle\n",
    "$$\n",
    "\n",
    "> **Warning!** Remember that every `Qobj` lives in a truncated Hilbert space. If the $\\alpha$ value is too large, the state will become a non-physical state because it will touch the high energy levels of the truncated Hilbert space\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "239e97ea",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense$$\\left(\\begin{array}{cc}0.726\\\\0.581\\\\0.329\\\\0.152\\\\0.061\\\\0.022\\\\0.007\\\\0.002\\\\6.048\\times10^{ -4 }\\\\1.713\\times10^{ -4 }\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense\n",
       "Qobj data =\n",
       "[[7.26149037e-01]\n",
       " [5.80919230e-01]\n",
       " [3.28617541e-01]\n",
       " [1.51781941e-01]\n",
       " [6.07127755e-02]\n",
       " [2.17212702e-02]\n",
       " [7.09408207e-03]\n",
       " [2.14540751e-03]\n",
       " [6.04780881e-04]\n",
       " [1.71316475e-04]]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "alpha = 0.8\n",
    "coherent_state = qutip.coherent(N, alpha)\n",
    "\n",
    "coherent_state"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "47b92030",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense$$\\left(\\begin{array}{cc}0.726\\\\0.581\\\\0.329\\\\0.152\\\\0.061\\\\0.022\\\\0.007\\\\0.002\\\\6.424\\times10^{ -4 }\\\\0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [1]], shape=(10, 1), type='ket', dtype=Dense\n",
       "Qobj data =\n",
       "[[7.26149037e-01]\n",
       " [5.80919230e-01]\n",
       " [3.28617541e-01]\n",
       " [1.51781939e-01]\n",
       " [6.07127959e-02]\n",
       " [2.17211016e-02]\n",
       " [7.09526842e-03]\n",
       " [2.13822331e-03]\n",
       " [6.42436781e-04]\n",
       " [0.00000000e+00]]"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "a * coherent_state / alpha"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "28f50106",
   "metadata": {},
   "source": [
    "### Spin States\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "e8b830b2",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[2], [1]], shape=(2, 1), type='ket', dtype=Dense$$\\left(\\begin{array}{cc}0\\\\1\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[2], [1]], shape=(2, 1), type='ket', dtype=Dense\n",
       "Qobj data =\n",
       "[[0.]\n",
       " [1.]]"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "qutip.spin_state(0.5, -1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "832c4dd7",
   "metadata": {},
   "source": [
    "### Density Matrices\n",
    "\n",
    "Quantum states can also be represented using density matrices, which is useful for describing mixed states.\n",
    "\n",
    "#### Creating a Density Matrix\n",
    "Let's convert our superposition state into a density matrix.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "44a015a7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Density matrix of the superposition state:\n"
     ]
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [10]], shape=(10, 10), type='oper', dtype=Dense, isherm=True$$\\left(\\begin{array}{cc}0.500 & 0.500 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0.500 & 0.500 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [10]], shape=(10, 10), type='oper', dtype=Dense, isherm=True\n",
       "Qobj data =\n",
       "[[0.5 0.5 0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.5 0.5 0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]\n",
       " [0.  0.  0.  0.  0.  0.  0.  0.  0.  0. ]]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Creating a density matrix from a state\n",
    "density_matrix = superposition_state * superposition_state.dag()  # Outer product\n",
    "\n",
    "print(\"Density matrix of the superposition state:\")\n",
    "density_matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0e8e2118",
   "metadata": {},
   "source": [
    "### Partial Trace\n",
    "\n",
    "The **partial trace** over a subsystem, say $B$, of a composite system $AB$, mathematically expresses as \"tracing out\" $B$, leaving the reduced state of $A$. For a bipartite state $\\rho_{AB}$, the partial trace over $B$ is:\n",
    "\n",
    "$$\n",
    "\\text{Tr}_B(\\hat{\\rho}_{AB}) = \\sum_{i \\in \\mathcal{H}_B} \\langle i| \\hat{\\rho}_{AB} |i\\rangle\n",
    "$$\n",
    "\n",
    "where $\\{|i\\rangle\\}$ forms a complete basis for subsystem $B$.\n",
    "\n",
    "Let's try it with an entangled Bell's state between two qubits:\n",
    "\n",
    "$$\n",
    "\\vert \\phi^+ \\rangle = \\frac{1}{\\sqrt{2}} \\left( \\vert 0, 0 \\rangle + \\vert 1, 1 \\rangle \\right)\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "53af30df",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=Dense, isherm=True$$\\left(\\begin{array}{cc}0.500 & 0\\\\0 & 0.500\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[2], [2]], shape=(2, 2), type='oper', dtype=Dense, isherm=True\n",
       "Qobj data =\n",
       "[[0.5 0. ]\n",
       " [0.  0.5]]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Bell state between two qubits\n",
    "phi_plus = ( qutip.tensor(qutip.spin_state(1/2, -1), qutip.spin_state(1/2, -1)) + qutip.tensor(qutip.spin_state(1/2, 1), qutip.spin_state(1/2, 1)) ).unit()\n",
    "\n",
    "# Reduced density matrix of the first qubit\n",
    "rho_1 = qutip.ptrace(phi_plus, 1)\n",
    "rho_1"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e718cdfe",
   "metadata": {},
   "source": [
    "We now apply the partial trace to a more complicated state, that is composed by two bosonic modes and two spins $\\vert j_1, m_1 \\rangle$ and $\\vert j_2, m_2 \\rangle$, with $j_1 = 1$ and $j_2 = \\frac{1}{2}$, $m_1 = 0$, and $m_2 = 1$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "b0640dad",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10, 10, 3], [10, 10, 3]], shape=(300, 300), type='oper', dtype=Dense, isherm=True$$\\left(\\begin{array}{cc}0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\\\vdots & \\vdots & \\vdots & \\vdots & \\vdots & \\ddots & \\vdots & \\vdots & \\vdots & \\vdots & \\vdots\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10, 10, 3], [10, 10, 3]], shape=(300, 300), type='oper', dtype=Dense, isherm=True\n",
       "Qobj data =\n",
       "[[0. 0. 0. ... 0. 0. 0.]\n",
       " [0. 0. 0. ... 0. 0. 0.]\n",
       " [0. 0. 0. ... 0. 0. 0.]\n",
       " ...\n",
       " [0. 0. 0. ... 0. 0. 0.]\n",
       " [0. 0. 0. ... 0. 0. 0.]\n",
       " [0. 0. 0. ... 0. 0. 0.]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10, 3], [10, 3]], shape=(30, 30), type='oper', dtype=Dense, isherm=True$$\\left(\\begin{array}{cc}0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 1 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\\\vdots & \\vdots & \\vdots & \\vdots & \\vdots & \\ddots & \\vdots & \\vdots & \\vdots & \\vdots & \\vdots\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & \\cdots & 0 & 0 & 0 & 0 & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10, 3], [10, 3]], shape=(30, 30), type='oper', dtype=Dense, isherm=True\n",
       "Qobj data =\n",
       "[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.\n",
       "  0. 0. 0. 0. 0. 0.]]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/latex": [
       "Quantum object: dims=[[10], [10]], shape=(10, 10), type='oper', dtype=Dense, isherm=True$$\\left(\\begin{array}{cc}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\end{array}\\right)$$"
      ],
      "text/plain": [
       "Quantum object: dims=[[10], [10]], shape=(10, 10), type='oper', dtype=Dense, isherm=True\n",
       "Qobj data =\n",
       "[[0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 1. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
       " [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]]"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "N = 10\n",
    "j1 = 1\n",
    "j2 = 1/2\n",
    "m1 = 0\n",
    "m2 = 1\n",
    "\n",
    "\n",
    "psi = qutip.tensor(qutip.fock(N, 3), qutip.fock(N, 1), qutip.spin_state(j1, 0), qutip.spin_state(j2, 1))\n",
    "\n",
    "# Trace only the second spin state\n",
    "rho_0 = qutip.ptrace(psi, [0, 1, 2])\n",
    "display(rho_0)\n",
    "\n",
    "# Trace only the first bosonic mode and the second spin state\n",
    "rho_1 = qutip.ptrace(psi, [1, 2])\n",
    "display(rho_1)\n",
    "\n",
    "# Trace all except the second bosonic mode\n",
    "rho_2 = qutip.ptrace(psi, [1])\n",
    "rho_2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "09a5a42e",
   "metadata": {},
   "source": [
    "## Eigenstates and Eigenvalues\n",
    "\n",
    "The eigenstates and eigenvalues of a system or an operator provide crucial insights into its properties. Let's explore how to calculate these in QuTiP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "393256b0",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Eigenvalues of Sigma Z:\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array([-1.,  1.])"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "Eigenstates of Sigma Z:\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "array([Quantum object: dims=[[2], [1]], shape=(2, 1), type='ket', dtype=Dense\n",
       "       Qobj data =\n",
       "       [[ 0.]\n",
       "        [-1.]]                                                               ,\n",
       "       Quantum object: dims=[[2], [1]], shape=(2, 1), type='ket', dtype=Dense\n",
       "       Qobj data =\n",
       "       [[-1.]\n",
       "        [-0.]]                                                               ],\n",
       "      dtype=object)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Example: Eigenstates and eigenvalues of Pauli Z\n",
    "eigenvalues, eigenstates = sigma_z.eigenstates()\n",
    "\n",
    "print(\"Eigenvalues of Sigma Z:\")\n",
    "display(eigenvalues)\n",
    "print(\"\\n\")\n",
    "print(\"Eigenstates of Sigma Z:\")\n",
    "display(eigenstates)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6fef8b86",
   "metadata": {},
   "source": [
    "\n",
    "## Computing Expectation Values\n",
    "\n",
    "The expectation value of an operator provides insight into the average outcome of a quantum measurement. For a quantum state $|\\psi\\rangle$ and an operator $\\hat{O}$, the expectation value is given by:\n",
    "\n",
    "$$\\langle \\hat{O} \\rangle = \\langle\\psi| \\hat{O} |\\psi\\rangle$$\n",
    "\n",
    "Expectation values are crucial for predicting measurable quantities in quantum mechanics. Let's compute the expectation value of the number operator $\\hat{n} = \\hat{a}^\\dagger \\hat{a}$ for a coherent state, which represents a quantum state closest to a classical harmonic oscillator.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "ffa66fa0",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Expectation value of the number operator for |psi>:\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "0.6399999989126255"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "The squared modulus of alpha is:\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "0.6400000000000001"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define the coherent state |psi> with alpha=2\n",
    "alpha = 0.8\n",
    "psi = qutip.coherent(N, alpha)\n",
    "\n",
    "# Define the number operator n = a.dag() * a\n",
    "n = a.dag() * a\n",
    "\n",
    "# Compute the expectation value of n for the state |psi>\n",
    "expectation_value_n = qutip.expect(n, psi)\n",
    "\n",
    "print(\"Expectation value of the number operator for |psi>:\")\n",
    "display(expectation_value_n)\n",
    "print(\"\\n\")\n",
    "print(\"The squared modulus of alpha is:\")\n",
    "display(abs(alpha) ** 2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e6cd7c28",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "physics",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}