From 87b964db3d1d745f9b1f50a06dbdd34aad6b8db9 Mon Sep 17 00:00:00 2001
From: Kees van Kempen <k.vankempen@student.science.ru.nl>
Date: Mon, 4 Apr 2022 20:40:20 +0200
Subject: [PATCH] Final: Do task 4.5

---
 ...l - Tight-binding propagation method.ipynb | 86 ++++++++++++++++---
 1 file changed, 73 insertions(+), 13 deletions(-)

diff --git a/Final/Final - Tight-binding propagation method.ipynb b/Final/Final - Tight-binding propagation method.ipynb
index 334d3e7..3918ebc 100644
--- a/Final/Final - Tight-binding propagation method.ipynb	
+++ b/Final/Final - Tight-binding propagation method.ipynb	
@@ -1693,7 +1693,7 @@
   },
   {
    "cell_type": "code",
-   "execution_count": null,
+   "execution_count": 20,
    "metadata": {
     "deletable": false,
     "nbgrader": {
@@ -1711,7 +1711,31 @@
    "outputs": [],
    "source": [
     "def getU_TZ(tau, H):\n",
-    "    # TODO: Add docstring.\n",
+    "    \"\"\"\n",
+    "    Calculates the time-propagation matrix for time-step tau using\n",
+    "    the Trotter-Suzuki algorithm for a known shape of the hamiltonian\n",
+    "    \n",
+    "    Args:\n",
+    "        tau: time-step numeric to calculate the time-propagation matrix at\n",
+    "        H:   n by n array representing the hamiltonian matrix\n",
+    "\n",
+    "    Returns:\n",
+    "        The time-propagation matrix U(tau).\n",
+    "    \n",
+    "    Note:\n",
+    "        Only the size of H is taken into account, as the decomposition is\n",
+    "        pre-defined.\n",
+    "    \"\"\"\n",
+    "    \n",
+    "    # We assume that H is n by n in the form with constant hopping parameter\n",
+    "    # t on the upper and lower diagonal, and the rest zero.\n",
+    "    n = len(H)\n",
+    "    assert len(H[0] == n)\n",
+    "    \n",
+    "    t = H[0][1]\n",
+    "    H_tb = (np.eye(n, n, -1) + np.eye(n, n, 1))*t\n",
+    "    assert np.allclose(H, H_tb)\n",
+    "    \n",
     "    \n",
     "    # First we calculate the exponent of the block.\n",
     "    exp_block = np.array([\n",
@@ -1721,23 +1745,59 @@
     "    \n",
     "    # Now we use the block to calculate the exponent of matrices\n",
     "    # H_1 and H_2.\n",
-    "    n = len(H)\n",
-    "    exp_H_1 = np.block([\n",
-    "        [1, 0],\n",
-    "        [0, np.kron(np.eye(n//2, dtype=int), exp_block)]\n",
-    "    ])\n",
-    "    exp_H_2 = np.block([\n",
-    "        [1,               0],\n",
-    "        [0, np.kron(np.eye(n//2, dtype=int), exp_block)]\n",
-    "    ])\n",
     "    \n",
-    "    # YOUR CODE HERE\n",
-    "    raise NotImplementedError()\n",
+    "    # If n is even.\n",
+    "    if n % 2 == 0:\n",
+    "        num = n//2\n",
+    "        exp_H_1 = np.kron(np.eye(num, dtype=int), exp_block)\n",
+    "        \n",
+    "        num = n//2 - 1\n",
+    "        exp_H_2 = np.block([\n",
+    "            [1,                    np.zeros((1, num*2)),                       0                   ],\n",
+    "            [np.zeros((num*2, 1)), np.kron(np.eye(num, dtype=int), exp_block), np.zeros((num*2, 1))],\n",
+    "            [0,                    np.zeros((1, num*2)),                       1                   ]\n",
+    "        ])\n",
+    "    \n",
+    "    # If n is odd.\n",
+    "    else:\n",
+    "        num = n//2\n",
+    "        exp_H_1 = np.block([\n",
+    "            [np.kron(np.eye(num, dtype=int), exp_block), np.zeros((num*2, 1))],\n",
+    "            [np.zeros((1, num*2)),                       1                   ]\n",
+    "        ])\n",
+    "\n",
+    "        num = n//2\n",
+    "        exp_H_2 = np.block([\n",
+    "            [1,                    np.zeros((1, num*2))                      ],\n",
+    "            [np.zeros((num*2, 1)), np.kron(np.eye(num, dtype=int), exp_block)]\n",
+    "        ])\n",
+    "\n",
+    "    U_TZ = exp_H_1@exp_H_2\n",
+    "    return U_TZ\n",
     "\n",
     "# Yann mentions again that this is slightly different wrong what\n",
     "# is in the slides/lecture."
    ]
   },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "n = 10\n",
+    "sigma = .25\n",
+    "tau = 1.5\n",
+    "H = TBHamiltonian(n, sigma)\n",
+    "U = getU_TZ(tau, H)\n",
+    "\n",
+    "exp_block = np.array([\n",
+    "        [    np.cos(tau*t),    -1j*np.sin(tau*t)],\n",
+    "        [-1j*np.sin(tau*t),        np.cos(tau*t)]\n",
+    "    ])\n",
+    "exp_block.shape"
+   ]
+  },
   {
    "cell_type": "markdown",
    "metadata": {