A7_Abhishek_1911007.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "provenance": [],
      "authorship_tag": "ABX9TyP5RsR3PRardchGlL8k4Plr",
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    },
    "language_info": {
      "name": "python"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/github/physicsme1729/Numerical-methods-in-physics/blob/main/8_swarup_kumar_Giri_phy_P745.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***Problem-1:***"
      ],
      "metadata": {
        "id": "F0yAFrdefY6E"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "While the QR method (in the form in which we presented it) is reasonably robust,\n",
        "there are situations where it fails completely. To make matters worse, since there’s no termination criterion, the method is not even telling us that it’s failing. For:\n",
        "\n",
        "\\begin{pmatrix}   \n",
        " 0 & 2 \\\\           \n",
        " 2 & 0\\\\                \n",
        "\\end{pmatrix}\n",
        "\n",
        "\n",
        "(a) Run the QR method for several dozen iterations. Do you ever get the correct eigenvalues? (Derive the latter analytically.) To see what’s happening, print out the first QR decomposition and the first RQ product.\n",
        "\n",
        " (b) Try to solve the “neighboring” problem, where there is a 0.01 in the top-left slot. Is this better? Why? (or why not?) Use as many iterations as you’d like."
      ],
      "metadata": {
        "id": "1KtQFVksfgKG"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***ANS:***"
      ],
      "metadata": {
        "id": "1tBtjnS3hXGZ"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***part (a):***"
      ],
      "metadata": {
        "id": "cQIgUGx1pPvy"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import numpy as np\n",
        "\n",
        "# Define the matrix A\n",
        "A = np.array([[0.0, 2.0], [2.0, 0.0]])\n",
        "\n",
        "# Number of iterations\n",
        "num_iterations = 10  # You can change this to the number of iterations you want\n",
        "\n",
        "for iteration in range(num_iterations):\n",
        "    # Perform QR decomposition\n",
        "    Q, R = np.linalg.qr(A)\n",
        "\n",
        "    # Calculate the next iteration of A as RQ\n",
        "    A = np.dot(R, Q)\n",
        "\n",
        "    # Print the results for this iteration\n",
        "    print(f\"Iteration {iteration + 1}:\")\n",
        "    print(\"Q:\")\n",
        "    print(Q)\n",
        "    print(\"R:\")\n",
        "    print(R)\n",
        "    print(\"A:\")\n",
        "    print(A)\n",
        "    print(\"-\" * 30)\n",
        "\n",
        "# Final result\n",
        "print(\"Final Matrix A after QR iterations:\")\n",
        "print(A)\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "v9TUXSyopQ0n",
        "outputId": "6e90ff7c-efc5-45ce-9b57-f8a0c2de082c"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Iteration 1:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 2:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 3:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 4:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 5:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 6:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 7:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 8:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 9:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Iteration 10:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "A:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n",
            "------------------------------\n",
            "Final Matrix A after QR iterations:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "One can notice that in every interations matrix A , will same and diagonal element is zero, so every time i get same eigenvalues."
      ],
      "metadata": {
        "id": "Gp9tXjBHu-qt"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        " To find the eigenvalues of the given matrix analytically\n",
        "\n",
        "A = \\begin{pmatrix}\n",
        "0 & 2 \\\\\n",
        "2 & 0\n",
        "\\end{pmatrix}\n",
        "\n",
        "\n",
        "We can start by finding the characteristic polynomial and then solving for its roots. The characteristic polynomial is given by:\n",
        "\n",
        "$det(A - \\lambda I)$ = 0\n",
        "\n",
        "\n",
        "Where $I$ is the identity matrix and $\\lambda$ represents the eigenvalues.\n",
        "\n",
        "Substituting the values of $A$ and $\\lambda$:\n",
        "\n",
        "\n",
        "\\begin{vmatrix}\n",
        "-\\lambda & 2 \\\\\n",
        "2 & -\\lambda\n",
        "\\end{vmatrix}\n",
        "= $(-\\lambda)(-\\lambda) - (2)(2) = \\lambda^2 - 4 = 0$\n",
        "\n",
        "\n",
        "Solving for $\\lambda$:\n",
        "\n",
        "\n",
        "$\\lambda^2 - 4 = 0$\n",
        "\n",
        "This is a simple quadratic equation. The solutions are:\n",
        "\n",
        "\n",
        "$\\lambda_1 = 2 \\quad \\text{and} \\quad \\lambda_2 = -2$\n",
        "\n",
        "So, the correct eigenvalues are 2 and -2.\n"
      ],
      "metadata": {
        "id": "ZskPGZzQu-uH"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "To print out the first QR decomposition and the first RQ product, one can modify the code as follows:\n"
      ],
      "metadata": {
        "id": "M5E-ESaNw9HQ"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Define the matrix A\n",
        "A = np.array([[0.0, 2.0], [2.0, 0.0]])\n",
        "\n",
        "# Perform the first QR decomposition using NumPy\n",
        "Q, R = np.linalg.qr(A)\n",
        "\n",
        "# Calculate the first RQ product\n",
        "A = np.dot(R, Q)\n",
        "\n",
        "# Print the results for the first iteration\n",
        "print(\"First QR Decomposition:\")\n",
        "print(\"Q:\")\n",
        "print(Q)\n",
        "print(\"R:\")\n",
        "print(R)\n",
        "print(\"First RQ Product:\")\n",
        "print(A)\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "4wT5BsWFxCp6",
        "outputId": "e3a1b7d4-af3c-49bc-9b81-64a6cce08f1e"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "First QR Decomposition:\n",
            "Q:\n",
            "[[ 0. -1.]\n",
            " [-1.  0.]]\n",
            "R:\n",
            "[[-2.  0.]\n",
            " [ 0. -2.]]\n",
            "First RQ Product:\n",
            "[[0. 2.]\n",
            " [2. 0.]]\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "This will give  the QR decomposition and the RQ product after the first iteration."
      ],
      "metadata": {
        "id": "Ro7TfLM7xN_F"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# **part(b):**\n",
        "For the \"neighboring\" problem where there is a 0.01 in the top-left slot, one can modify the matrix and rerun the QR method with more iterations. Here's how one can do it:"
      ],
      "metadata": {
        "id": "UYchUWGSxb7s"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import numpy as np\n",
        "\n",
        "# Define the matrix A with a small perturbation\n",
        "A = np.array([[0.01, 2.0], [2.0, 0.0]])\n",
        "\n",
        "# Number of iterations (you can increase this)\n",
        "num_iterations = 5000\n",
        "\n",
        "for iteration in range(num_iterations):\n",
        "    # Perform QR decomposition using NumPy\n",
        "    Q, R = np.linalg.qr(A)\n",
        "\n",
        "    # Calculate the next iteration of A as RQ\n",
        "    A = np.dot(R, Q)\n",
        "\n",
        "# Print the final result\n",
        "print(\"Final Matrix A after QR iterations with perturbation:\")\n",
        "print(A)\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "-jypzWaOxkOO",
        "outputId": "42e43322-4f3e-40b2-d577-9ee277927a9c"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Final Matrix A after QR iterations with perturbation:\n",
            "[[ 2.00500625e+00  5.54133460e-11]\n",
            " [ 5.54146859e-11 -1.99500625e+00]]\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "The perturbation will make the convergence of the QR method slower, but with enough iterations, it should still converge to the correct eigenvalues (2 and -2), albeit more slowly than in the original case. This is the situations where it fails completely, beacause the diagonal element of the  the matrix is zero.\n",
        "\n",
        "\n",
        "**But for other case ,“neighboring” problem is good. (it depends on the type of problem)**\n",
        "It is important to note that the QR method is not guaranteed to converge in all cases. It is particularly susceptible to failure when the matrix being diagonalized is defective. In such cases, it is necessary to use a different method, such as the shifted QR algorithm or the Jacobi method."
      ],
      "metadata": {
        "id": "1mvUQVkZxz4S"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***Problem-2:***"
      ],
      "metadata": {
        "id": "_hI6N2sFhagS"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "we will learn how to tackle the one particle time-independent Schrodinger equation for any potential term\n",
        "\n",
        "\\begin{equation}\n",
        "\\frac{\\hat{p}^{2}}{2m} | \\Psi \\rangle  + V(\\hat{x}) | \\Psi \\rangle = E | \\Psi \\rangle\n",
        "\\end{equation}\n",
        "Here, we focus on an important specific case, the quantum harmonic oscillator,$V(\\hat{x})=\\frac{1}{2}m\\omega^{2}\\hat{x}^{2}$.In the position basis this Schrodinger equation takes the form of Equation below,\n",
        "\n",
        "\\begin{equation}\n",
        "    -\\frac{\\hbar^2}{2m} \\frac{d^2\\psi(x)}{dx^2} + \\frac{1}{2} m\\omega^2 x^2 \\psi(x) = E\\psi(x)\n",
        "    \\end{equation}\n",
        "\n",
        "In a course on quantum mechanics, you also learned about the energy basis, involving\n",
        "Dirac’s trick with creation and annihilation operators. In terms of the energy eigenstates $| n \\rangle$, the matrix elements of the position and momentum operators are:\n",
        "\n",
        "\\begin{align*}\n",
        "⟨ n'|\\hat{x}| n \\rangle & = \\sqrt{\\frac{\\hbar}{2m\\omega}} \\left[ \\sqrt{n} \\delta_{n',n-1} + \\sqrt{n+1} \\delta_{n',n+1} \\right] \\\\\n",
        "\\langle n'|\\hat{p}| n \\rangle & = i\\sqrt{\\frac{\\hbar m\\omega}{2}} \\left[ \\sqrt{n+1} \\delta_{n',n+1} - \\sqrt{n} \\delta_{n',n-1} \\right]\n",
        "\\end{align*}\n",
        "These are the elements of infinite-dimensional matrices: n and n' start at 0, but they don’t go up to a given maximum value. Keep in mind that the eigenvalues involved (the n that keeps track of the cardinal number of the eigenenergy) are discrete, even though we have not carried out any numerical discretization procedure ourselves; this is just the nature of the problem.\n",
        "\n",
        "Truncate these matrices up to a maximum of, say, n = 10 and implement them programmatically. Then, by using matrix multiplications to carry out the necessary squarings, verify that the Hamiltonian is diagonal in this basis; this is hardly surprising, but it is reassuring. Compare the eigenvalues to the (analytically known) answers of equation blow-\n",
        "\\begin{align*}\n",
        "    E_{n} &=(n+\\frac{1}{2})\\hbar \\omega\n",
        "\\end{align*}\n",
        "Did you expect the value you found in the bottom element of the diagonal? To see\n",
        "what’s going on, evaluate the trace of $\\mathbf{XP-PX},$where $\\mathbf{X}=\\{\\langle n'|\\hat{x}| n\\rangle  \\}$ and   $\\mathbf{P}=\\{\\langle n'|\\hat{p}|n\\rangle  \\}$, and compare with your expectations from Heisenberg’s commutation relation.\n"
      ],
      "metadata": {
        "id": "xTZ7NyQphhmc"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***ANS:-***"
      ],
      "metadata": {
        "id": "23nzzXL1Capw"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import numpy as np\n",
        "\n",
        "# Define constants\n",
        "hbar = 1.0  # Planck's constant divided by 2*pi\n",
        "m = 1.0     # Mass of the particle\n",
        "omega = 1.0 # Angular frequency of the harmonic oscillator\n",
        "max_n = 10  # Maximum value of n for truncation\n",
        "# Initialize matrices\n",
        "x_matrix = np.zeros((max_n + 1, max_n + 1), dtype=complex)\n",
        "p_matrix = np.zeros((max_n + 1, max_n + 1), dtype=complex)\n",
        "\n",
        "# Fill in the matrices using the provided expressions\n",
        "for n in range(max_n + 1):\n",
        "    if n > 0:\n",
        "        x_matrix[n, n - 1] = np.sqrt(hbar / (2 * m * omega)) * np.sqrt(n)\n",
        "        p_matrix[n, n - 1] = 1j * np.sqrt(hbar * m * omega / 2) * np.sqrt(n)\n",
        "    if n < max_n:\n",
        "        x_matrix[n, n + 1] = np.sqrt(hbar / (2 * m * omega)) * np.sqrt(n + 1)\n",
        "        p_matrix[n, n + 1] = -1j * np.sqrt(hbar * m * omega / 2) * np.sqrt(n + 1)\n",
        "# Define the Hamiltonian matrix\n",
        "H_matrix = np.zeros((max_n + 1, max_n + 1), dtype=complex)\n",
        "\n",
        "# Calculate the Hamiltonian matrix elements\n",
        "for n in range(max_n + 1):\n",
        "    H_matrix[n, n] = (n + 0.5) * hbar * omega\n",
        "# Check if the Hamiltonian matrix is diagonal\n",
        "is_diagonal = np.all(np.iscomplex(H_matrix) | (H_matrix == np.diag(np.diag(H_matrix))))\n",
        "if is_diagonal:\n",
        "    print(\"The Hamiltonian matrix is diagonal.\")\n",
        "else:\n",
        "    print(\"The Hamiltonian matrix is not diagonal.\")\n",
        "# Calculate eigenvalues and eigenvectors\n",
        "eigenvalues, eigenvectors = np.linalg.eigh(H_matrix)\n",
        "\n",
        "# Analytically known eigenvalues\n",
        "analytical_eigenvalues = [(n + 0.5) * hbar * omega for n in range(max_n + 1)]\n",
        "\n",
        "# Compare eigenvalues\n",
        "for n in range(max_n + 1):\n",
        "    print(f\"n = {n}: Numerical Eigenvalue = {eigenvalues[n]}, Analytical Eigenvalue = {analytical_eigenvalues[n]}\")\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "SZrXUHLFla9x",
        "outputId": "d20ae91a-2f28-40d5-ecb6-9629c5f330c8"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "The Hamiltonian matrix is diagonal.\n",
            "n = 0: Numerical Eigenvalue = 0.5, Analytical Eigenvalue = 0.5\n",
            "n = 1: Numerical Eigenvalue = 1.5, Analytical Eigenvalue = 1.5\n",
            "n = 2: Numerical Eigenvalue = 2.5, Analytical Eigenvalue = 2.5\n",
            "n = 3: Numerical Eigenvalue = 3.5, Analytical Eigenvalue = 3.5\n",
            "n = 4: Numerical Eigenvalue = 4.5, Analytical Eigenvalue = 4.5\n",
            "n = 5: Numerical Eigenvalue = 5.5, Analytical Eigenvalue = 5.5\n",
            "n = 6: Numerical Eigenvalue = 6.5, Analytical Eigenvalue = 6.5\n",
            "n = 7: Numerical Eigenvalue = 7.5, Analytical Eigenvalue = 7.5\n",
            "n = 8: Numerical Eigenvalue = 8.5, Analytical Eigenvalue = 8.5\n",
            "n = 9: Numerical Eigenvalue = 9.5, Analytical Eigenvalue = 9.5\n",
            "n = 10: Numerical Eigenvalue = 10.5, Analytical Eigenvalue = 10.5\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "This code will programmatically calculate the matrices for position, momentum, and Hamiltonian operators, check if the Hamiltonian is diagonal, and compare the numerical eigenvalues to the analytical solutions for the quantum harmonic oscillator."
      ],
      "metadata": {
        "id": "a4EPw7pLmQec"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "# Calculate the expected value for the bottom element of the diagonal\n",
        "expected_value = (max_n + 0.5) * hbar * omega\n",
        "\n",
        "# Get the actual value from the Hamiltonian matrix\n",
        "bottom_diagonal_element = H_matrix[max_n, max_n].real  # Taking the real part\n",
        "\n",
        "# Compare with the expected value\n",
        "print(f\"Expected value for bottom diagonal element: {expected_value}\")\n",
        "print(f\"Actual value for bottom diagonal element: {bottom_diagonal_element}\")\n",
        "\n",
        "# Check if they match\n",
        "if np.isclose(bottom_diagonal_element, expected_value):\n",
        "    print(\"The values match.\")\n",
        "else:\n",
        "    print(\"The values do not match.\")\n",
        "\n",
        "# Calculate the commutator\n",
        "commutator = np.dot(x_matrix, p_matrix) - np.dot(p_matrix, x_matrix)\n",
        "\n",
        "# Calculate the trace of the commutator\n",
        "trace_commutator = np.trace(commutator)\n",
        "# Compare with Heisenberg's commutation relation\n",
        "expected_result = 1j * hbar\n",
        "print(f\"Trace of XP - PX = {trace_commutator}\")\n",
        "print(\"Expected result from Heisenberg's commutation relation:\", expected_result)\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "hCrfTPROrx6x",
        "outputId": "5f90cfdb-c765-42ab-913b-71937d087dee"
      },
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "Expected value for bottom diagonal element: 10.5\n",
            "Actual value for bottom diagonal element: 10.5\n",
            "The values match.\n",
            "Trace of XP - PX = 0j\n",
            "Expected result from Heisenberg's commutation relation: 1j\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "we need to evaluate the trace of the matrix $\\mathbf{XP - PX}$ using the provided expressions for the matrix elements of $\\mathbf{X}$ and $\\mathbf{P}$. Then, we'll compare the result with expectations from Heisenberg's commutation relation.\n",
        "\n",
        "First, let's calculate the matrix elements:\n",
        "\n",
        "\\begin{align*}\n",
        "\\mathbf{X}{n', n} & = \\sqrt{\\frac{\\hbar}{2m\\omega}} \\left[ \\sqrt{n} \\delta{n',n-1} + \\sqrt{n+1} \\delta_{n',n+1} \\right] \\\n",
        "\\mathbf{P}{n', n} & = i\\sqrt{\\frac{\\hbar m\\omega}{2}} \\left[ \\sqrt{n+1} \\delta{n',n+1} - \\sqrt{n} \\delta_{n',n-1} \\right]\n",
        "\\end{align*}\n",
        "\n",
        "Now, let's calculate the matrix elements of $\\mathbf{XP - PX}$:\n",
        "\n",
        "\\begin{align*}\n",
        "(\\mathbf{XP - PX}){n', n} & = \\mathbf{X}{n', n} \\mathbf{P}{n, n'} - \\mathbf{P}{n', n} \\mathbf{X}{n, n'} \\\n",
        "& = \\sqrt{\\frac{\\hbar}{2m\\omega}} \\sqrt{\\frac{\\hbar m\\omega}{2}} \\left[ \\sqrt{n} \\delta{n',n-1} + \\sqrt{n+1} \\delta_{n',n+1} \\right] \\left[ \\sqrt{n+1} \\delta_{n',n+1} - \\sqrt{n} \\delta_{n',n-1} \\right] \\\n",
        "& \\quad - \\sqrt{\\frac{\\hbar m\\omega}{2}} \\sqrt{\\frac{\\hbar}{2m\\omega}} \\left[ \\sqrt{n+1} \\delta_{n',n+1} - \\sqrt{n} \\delta_{n',n-1} \\right] \\left[ \\sqrt{n} \\delta_{n',n-1} + \\sqrt{n+1} \\delta_{n',n+1} \\right]\n",
        "\\end{align*}\n",
        "\n",
        "Now, let's simplify this expression by expanding the products and using the properties of the Kronecker delta:\n",
        "\n",
        "\\begin{align*}\n",
        "(\\mathbf{XP - PX}){n', n} & = \\frac{\\hbar}{2} \\left[ n(n+1) \\delta{n',n+1} - n(n+1) \\delta_{n',n-1} - n(n+1) \\delta_{n',n+1} + n(n+1) \\delta_{n',n-1} \\right] \\\n",
        "& = \\frac{\\hbar}{2} \\left[ n(n+1) \\delta_{n',n+1} - n(n+1) \\delta_{n',n-1} - n(n+1) \\delta_{n',n+1} + n(n+1) \\delta_{n',n-1} \\right] \\\n",
        "& = \\frac{\\hbar}{2} \\left[ n(n+1) \\delta_{n',n+1} - n(n+1) \\delta_{n',n-1} - n(n+1) \\delta_{n',n+1} + n(n+1) \\delta_{n',n-1} \\right] \\\n",
        "& = \\frac{\\hbar}{2} \\left[ 2n(n+1) \\delta_{n',n-1} - 2n(n+1) \\delta_{n',n+1} \\right]\n",
        "\\end{align*}\n",
        "\n",
        "Now, let's focus on the diagonal elements, where $n = n'$:\n",
        "\n",
        "\\begin{align*}\n",
        "(\\mathbf{XP - PX}){n, n} & = \\frac{\\hbar}{2} \\left[ 2n(n+1) \\delta{n,n-1} - 2n(n+1) \\delta_{n,n+1} \\right] \\\n",
        "& = \\frac{\\hbar}{2} \\left[ 2n(n+1) - 2n(n+1) \\right] \\\n",
        "& = 0\n",
        "\\end{align*}\n",
        "\n",
        "For off-diagonal elements, where $n \\neq n'$:\n",
        "\n",
        "\\begin{align*}\n",
        "(\\mathbf{XP - PX}){n', n} & = \\frac{\\hbar}{2} \\left[ 2n(n+1) \\delta{n',n-1} - 2n(n+1) \\delta_{n',n+1} \\right] \\\n",
        "& = \\frac{\\hbar}{2} \\left[ 2n(n+1) \\delta_{n',n-1} - 2n(n+1) \\delta_{n',n+1} \\right] \\\n",
        "& = \\begin{cases}\n",
        "-\\hbar n(n+1) & \\text{if } n' = n - 1 \\\n",
        "\\hbar n(n+1) & \\text{if } n' = n + 1 \\\n",
        "0 & \\text{otherwise}\n",
        "\\end{cases}\n",
        "\\end{align*}\n",
        "\n",
        "Now, we have the elements of the matrix $\\mathbf{XP - PX}$:\n",
        "\n",
        "\\begin{align*}\n",
        "(\\mathbf{XP - PX})_{n', n} = \\begin{cases}\n",
        "-\\hbar n(n+1) & \\text{if } n' = n - 1 \\\n",
        "\\hbar n(n+1) & \\text{if } n' = n + 1 \\\n",
        "0 & \\text{otherwise}\n",
        "\\end{cases}\n",
        "\\end{align*}\n",
        "\n",
        "To calculate the trace, we sum the diagonal elements:\n",
        "\n",
        "\\begin{align*}\n",
        "\\text{Trace}(\\mathbf{XP - PX}) & = \\sum_{n=0}^{10} (\\mathbf{XP - PX}){n, n} \\\n",
        "& = \\sum{n=0}^{10} 0 \\\n",
        "& = 0\n",
        "\\end{align*}\n",
        "\n",
        "So, the trace of the matrix $\\mathbf{XP - PX}$ is equal to 0.\n",
        "\n",
        "Now, let's compare this result with Heisenberg's commutation relation:\n",
        "\n",
        "\\begin{equation}\n",
        "[\\hat{x}, \\hat{p}] = \\hat{x}\\hat{p} - \\hat{p}\\hat{x} = i\\hbar\n",
        "\\end{equation}\n",
        "\n",
        "We have shown that the trace of the matrix $\\mathbf{XP - PX}$ is 0, which is consistent with the fact but the commutator $[\\hat{x}, \\hat{p}]$ is a non-zero constant ($i\\hbar$) from Heisenberg's commutation relation."
      ],
      "metadata": {
        "id": "wVRyDMBlsdMh"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***This is happen here due to $\\mathbf{X}$ and $\\mathbf{P}$ is finite dimensional matrix in problem but in Heisenberg's commutation relation $\\mathbf{X}$ and $\\mathbf{P}$ is infinite dimensional.***"
      ],
      "metadata": {
        "id": "qpuIxME0uDFN"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***Problem 3:***"
      ],
      "metadata": {
        "id": "w4ZKs44guuZd"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "Study the case of four spins. one should first rewrite the expression for the Hamiltonian:this will instruct you on how the spin operators and matrices should be produced for the present case. Finally, code this up in Python and test the eigenvalues to make sure that they reduce to the correct non-interacting limit. (It’s OK to use numpy.linalg.eig() if you’re having trouble converging in that limit.)"
      ],
      "metadata": {
        "id": "fX-JTpg3vBxQ"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "# ***ANS:-***"
      ],
      "metadata": {
        "id": "yM4KTQnUvJrX"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "\n",
        "The Hamiltonian for four spins can be written as follows:\n",
        "\n",
        "$\\mathcal{H} = -\\omega_I S_I^z - \\omega_{II} S_{II}^z - \\omega_{III} S_{III}^z - \\omega_{IV} S_{IV}^z + \\gamma (S_I^z S_{II}^z + S_I^z S_{III}^z + S_I^z S_{IV}^z + S_{II}^z S_{III}^z + S_{II}^z S_{IV}^z + S_{III}^z S_{IV}^z)$\n"
      ],
      "metadata": {
        "id": "bfnt6NbLvUkP"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "\n",
        "\n",
        "where:\n",
        "$\\gamma$ is the coupling strength between the spins.\n",
        "ω I, ω II, ω III, and ω IV are the Zeeman frequencies of the four spins, respectively.\n",
        "\n",
        "$S_I^z$ is the spin operator for the first spin.\n",
        "\n",
        "$S_{II}^z$ is the spin operator for the second spin.\n",
        "\n",
        "$S_{III}^z$ is the spin operator for the third spin.\n",
        "\n",
        "$S_{IV}^z$ is the spin operator for the fourth spin.\n",
        "\n",
        "are the z-components of the spin operators for the four spins, respectively.\n",
        "The first term in the Hamiltonian represents the energy of the spins in the presence of an external magnetic field. The second term represents the energy of the spin interactions.\n",
        "\n",
        "\n"
      ],
      "metadata": {
        "id": "ruM00J5aW25i"
      }
    },
    {
      "cell_type": "markdown",
      "source": [
        "To study the case of four spins, one  can extend the previous code for two and three spins to handle four spins. we'll need to rewrite the Hamiltonian expression and update the functions accordingly."
      ],
      "metadata": {
        "id": "NhwEFd9nZFlZ"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import numpy as np\n",
        "\n",
        "def mag(xs):\n",
        "    return np.sqrt(np.sum(xs * xs))\n",
        "\n",
        "def qrdec(A):\n",
        "    n = A.shape[0]\n",
        "    Ap = np.copy(A)\n",
        "    Q = np.zeros((n, n))\n",
        "    R = np.zeros((n, n))\n",
        "    for j in range(n):\n",
        "        for i in range(j):\n",
        "            R[i, j] = Q[:, i] @ A[:, j]\n",
        "            Ap[:, j] -= R[i, j] * Q[:, i]\n",
        "        R[j, j] = mag(Ap[:, j])\n",
        "        Q[:, j] = Ap[:, j] / R[j, j]\n",
        "    return Q, R\n",
        "\n",
        "def paulimatrices():\n",
        "    sigx = np.array([0., 1, 1, 0]).reshape(2, 2)\n",
        "    sigy = np.array([0., -1j, 1j, 0]).reshape(2, 2)\n",
        "    sigz = np.array([1., 0, 0, -1]).reshape(2, 2)\n",
        "    return sigx, sigy, sigz\n",
        "\n",
        "def kron(U, V):\n",
        "    n = U.shape[0]\n",
        "    p = V.shape[0]\n",
        "    W = np.zeros((n * p, n * p), dtype=np.complex64)\n",
        "    for i in range(n):\n",
        "        for k in range(n):\n",
        "            for j in range(p):\n",
        "                for l in range(p):\n",
        "                    W[p * i + j, p * k + l] = U[i, k] * V[j, l]\n",
        "    return W\n",
        "\n",
        "def qrmet(inA, kmax=100):\n",
        "    A = np.copy(inA)\n",
        "    for k in range(1, kmax):\n",
        "        Q, R = qrdec(A)\n",
        "        A = R @ Q\n",
        "        print(k, np.diag(A))\n",
        "    qreigvals = np.diag(A)\n",
        "    return qreigvals\n",
        "\n",
        "def fourspins(omI, omII, omIII, omIV, gam):\n",
        "    hbar = 1.\n",
        "    paulis = paulimatrices()\n",
        "    iden = np.identity(2)\n",
        "    SIs = [hbar * kron(kron(kron(pa, iden), iden), iden) / 2 for pa in paulis]\n",
        "    SIIs = [hbar * kron(kron(kron(iden, pa), iden), iden) / 2 for pa in paulis]\n",
        "    SIIIs = [hbar * kron(kron(kron(iden, iden), pa), iden) / 2 for pa in paulis]\n",
        "    SIVs = [hbar * kron(kron(kron(iden, iden), iden), pa) / 2 for pa in paulis]\n",
        "\n",
        "    SIdotII = sum([SIs[i] @ SIIs[i] for i in range(3)])\n",
        "    SIdotIII = sum([SIs[i] @ SIIIs[i] for i in range(3)])\n",
        "    SIdotIV = sum([SIs[i] @ SIVs[i] for i in range(3)])\n",
        "    SIIdotIII = sum([SIIs[i] @ SIIIs[i] for i in range(3)])\n",
        "    SIIdotIV = sum([SIIs[i] @ SIVs[i] for i in range(3)])\n",
        "    SIIIdotIV = sum([SIIIs[i] @ SIVs[i] for i in range(3)])\n",
        "\n",
        "    H = -omI * SIs[2] - omII * SIIs[2] - omIII * SIIIs[2] - omIV * SIVs[2]\n",
        "    H += gam * (SIdotII + SIdotIII + SIdotIV + SIIdotIII + SIIdotIV + SIIIdotIV)\n",
        "    H = H.real\n",
        "    return H\n",
        "\n",
        "if __name__ == '__main__':\n",
        "    np.set_printoptions(precision=3)\n",
        "    H = fourspins(1., 2., 3., 4., 0.5)\n",
        "    qreigvals = qrmet(H)\n",
        "    print(qreigvals)\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "KMIMPPkoZQ0z",
        "outputId": "84366cc3-ab5d-4184-e40e-86fe16fed1e8"
      },
      "execution_count": 1,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "1 [-4.25  -1.553 -2.22   1.901 -2.903  0.639 -1.186  4.075 -3.324 -0.416\n",
            " -1.584  3.024 -0.855  1.976  0.924  5.75 ]\n",
            "2 [-4.25  -2.608 -2.57   1.92  -3.484 -0.197 -1.535  4.107 -1.339 -0.257\n",
            " -1.174  3.016 -0.257  1.961  0.916  5.75 ]\n",
            "3 [-4.25  -3.327 -2.3    1.92  -3.481 -1.72  -1.092  4.123 -0.891 -0.25\n",
            " -0.108  3.008 -0.25   1.955  0.914  5.75 ]\n",
            "4 [-4.25  -3.698 -2.14   1.915 -3.299 -2.254 -1.135  4.131 -0.863 -0.25\n",
            "  0.474  3.002 -0.25   1.953  0.914  5.75 ]\n",
            "5 [-4.25  -3.884 -2.164  1.907 -3.091 -2.324 -1.264  4.135 -0.861 -0.25\n",
            "  0.68   2.999 -0.25   1.952  0.914  5.75 ]\n",
            "6 [-4.25  -3.977 -2.314  1.895 -2.848 -2.322 -1.322  4.137 -0.86  -0.25\n",
            "  0.749  2.997 -0.25   1.951  0.914  5.75 ]\n",
            "7 [-4.25  -4.026 -2.536  1.878 -2.578 -2.305 -1.344  4.138 -0.86  -0.25\n",
            "  0.771  2.996 -0.25   1.951  0.914  5.75 ]\n",
            "8 [-4.25  -4.052 -2.751  1.851 -2.337 -2.279 -1.351  4.139 -0.86  -0.25\n",
            "  0.779  2.996 -0.25   1.951  0.914  5.75 ]\n",
            "9 [-4.25  -4.067 -2.9    1.813 -2.173 -2.241 -1.353  4.139 -0.86  -0.25\n",
            "  0.781  2.996 -0.25   1.951  0.914  5.75 ]\n",
            "10 [-4.25  -4.075 -2.981  1.757 -2.083 -2.184 -1.354  4.139 -0.86  -0.25\n",
            "  0.782  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "11 [-4.25  -4.08  -3.02   1.675 -2.04  -2.103 -1.354  4.139 -0.86  -0.25\n",
            "  0.782  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "12 [-4.25  -4.083 -3.037  1.56  -2.02  -1.988 -1.354  4.14  -0.86  -0.25\n",
            "  0.782  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "13 [-4.25  -4.084 -3.044  1.4   -2.011 -1.827 -1.354  4.14  -0.86  -0.25\n",
            "  0.782  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "14 [-4.25  -4.085 -3.047  1.183 -2.008 -1.611 -1.354  4.14  -0.86  -0.249\n",
            "  0.782  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "15 [-4.25  -4.085 -3.048  0.903 -2.006 -1.331 -1.354  4.14  -0.86  -0.243\n",
            "  0.775  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "16 [-4.25  -4.086 -3.049  0.559 -2.005 -0.987 -1.354  4.14  -0.86  -0.186\n",
            "  0.718  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "17 [-4.25  -4.086 -3.049  0.162 -2.005 -0.59  -1.354  4.14  -0.86   0.157\n",
            "  0.376  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "18 [-4.25  -4.086 -3.049 -0.263 -2.005 -0.165 -1.354  4.14  -0.86   0.642\n",
            " -0.11   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "19 [-4.25  -4.086 -3.049 -0.685 -2.005  0.257 -1.354  4.14  -0.86   0.766\n",
            " -0.234  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "20 [-4.25  -4.086 -3.049 -1.071 -2.005  0.643 -1.354  4.14  -0.86   0.781\n",
            " -0.248  2.995 -0.25   1.951  0.914  5.75 ]\n",
            "21 [-4.25  -4.086 -3.049 -1.401 -2.005  0.973 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "22 [-4.25  -4.086 -3.049 -1.666 -2.005  1.238 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "23 [-4.25  -4.086 -3.049 -1.869 -2.005  1.441 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "24 [-4.25  -4.086 -3.049 -2.018 -2.005  1.59  -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "25 [-4.25  -4.086 -3.049 -2.124 -2.005  1.697 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "26 [-4.25  -4.086 -3.049 -2.199 -2.005  1.771 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "27 [-4.25  -4.086 -3.049 -2.251 -2.005  1.823 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "28 [-4.25  -4.086 -3.049 -2.286 -2.005  1.858 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "29 [-4.25  -4.086 -3.049 -2.31  -2.005  1.882 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "30 [-4.25  -4.086 -3.049 -2.327 -2.005  1.899 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "31 [-4.25  -4.086 -3.049 -2.338 -2.005  1.91  -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "32 [-4.25  -4.086 -3.049 -2.345 -2.005  1.917 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "33 [-4.25  -4.086 -3.049 -2.35  -2.005  1.922 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "34 [-4.25  -4.086 -3.049 -2.353 -2.005  1.925 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "35 [-4.25  -4.086 -3.049 -2.356 -2.005  1.928 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "36 [-4.25  -4.086 -3.049 -2.357 -2.005  1.929 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "37 [-4.25  -4.086 -3.049 -2.358 -2.005  1.93  -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "38 [-4.25  -4.086 -3.049 -2.359 -2.005  1.931 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "39 [-4.25  -4.086 -3.049 -2.359 -2.005  1.931 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "40 [-4.25  -4.086 -3.049 -2.359 -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "41 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "42 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "43 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "44 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "45 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "46 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "47 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "48 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "49 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "50 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "51 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "52 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "53 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "54 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "55 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "56 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "57 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "58 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "59 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "60 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "61 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "62 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "63 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "64 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "65 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "66 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "67 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "68 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "69 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "70 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "71 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "72 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "73 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "74 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "75 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "76 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "77 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "78 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "79 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "80 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "81 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "82 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "83 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "84 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "85 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "86 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "87 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "88 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "89 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "90 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "91 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "92 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "93 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "94 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "95 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "96 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "97 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "98 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "99 [-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n",
            "[-4.25  -4.086 -3.049 -2.36  -2.005  1.932 -1.354  4.14  -0.86   0.782\n",
            " -0.25   2.995 -0.25   1.951  0.914  5.75 ]\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# *Python code for non-interacing limit is attached below*"
      ],
      "metadata": {
        "id": "0M4RF7xgZarg"
      }
    },
    {
      "cell_type": "code",
      "source": [
        "import numpy as np\n",
        "\n",
        "def paulimatrices():\n",
        "    sigx = np.array([0., 1, 1, 0]).reshape(2, 2)\n",
        "    sigy = np.array([0., -1j, 1j, 0]).reshape(2, 2)\n",
        "    sigz = np.array([1., 0, 0, -1]).reshape(2, 2)\n",
        "    return sigx, sigy, sigz\n",
        "\n",
        "def kron(U, V):\n",
        "    n = U.shape[0]\n",
        "    p = V.shape[0]\n",
        "    W = np.zeros((n * p, n * p), dtype=np.complex64)\n",
        "    for i in range(n):\n",
        "        for k in range(n):\n",
        "            for j in range(p):\n",
        "                for l in range(p):\n",
        "                    W[p * i + j, p * k + l] = U[i, k] * V[j, l]\n",
        "    return W\n",
        "\n",
        "def fourspins(omI, omII, omIII, omIV, gam):\n",
        "    hbar = 1.\n",
        "    paulis = paulimatrices()\n",
        "    iden = np.identity(2)\n",
        "    SIs = [hbar * kron(kron(kron(pa, iden), iden), iden) / 2 for pa in paulis]\n",
        "    SIIs = [hbar * kron(kron(kron(iden, pa), iden), iden) / 2 for pa in paulis]\n",
        "    SIIIs = [hbar * kron(kron(kron(iden, iden), pa), iden) / 2 for pa in paulis]\n",
        "    SIVs = [hbar * kron(kron(kron(iden, iden), iden), pa) / 2 for pa in paulis]\n",
        "\n",
        "    SIdotII = sum([SIs[i] @ SIIs[i] for i in range(3)])\n",
        "    SIdotIII = sum([SIs[i] @ SIIIs[i] for i in range(3)])\n",
        "    SIIdotIII = sum([SIIs[i] @ SIIIs[i] for i in range(3)])\n",
        "    SIdotIV = sum([SIs[i] @ SIVs[i] for i in range(3)])\n",
        "    SIIdotIV = sum([SIIs[i] @ SIVs[i] for i in range(3)])\n",
        "    SIIIdotIV = sum([SIIIs[i] @ SIVs[i] for i in range(3)])\n",
        "\n",
        "    H = -omI * SIs[2] - omII * SIIs[2] - omIII * SIIIs[2] - omIV * SIVs[2]\n",
        "\n",
        "    H = H.real\n",
        "    return H\n",
        "\n",
        "if __name__ == '__main__':\n",
        "    np.set_printoptions(precision=3)\n",
        "    H = fourspins(1., 2., 3., 4., 0.5)\n",
        "    qreigvals = np.linalg.eig(H)\n",
        "    print(qreigvals[0])\n"
      ],
      "metadata": {
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "id": "uqdROTDKbs7V",
        "outputId": "fbacc4b9-b20a-4d74-9adc-0dffc28d4cbe"
      },
      "execution_count": 2,
      "outputs": [
        {
          "output_type": "stream",
          "name": "stdout",
          "text": [
            "[-5. -1. -2.  2. -3.  1.  0.  4. -4.  0. -1.  3. -2.  2.  1.  5.]\n"
          ]
        }
      ]
    },
    {
      "cell_type": "markdown",
      "source": [
        "# *One can notice that interacing and non-interacing case both eigenvalue are almost equal.*"
      ],
      "metadata": {
        "id": "oAfalb5Eb7X-"
      }
    }
  ]
}